4.2. Blocks and Other Pyomo Best Practices#
4.2.1. Learning Objectives#
Practice using blocks to organize Pyomo models with a hierachical structure (e.g., two-stage stochastic program)
Learn best practices for organizing Pyomo models and data
4.2.2. Blocks in Pyomo#
Reference: Chapter 8: Stuctured Modeling with Blocks, Pyomo – Optimization Modeling in Python, Bynum et al. (2021).
Blocks provide a convenient way to express hierachically-structured models in Pyomo. As an example, consider the special structure in the Farmer’s example for Stochastic Programming.
Asking ChatGPT to rewrite our model from Stochastic Programming with blocks gives the following code (which was then slightly modified):
from pyomo.environ import ConcreteModel, Var, Objective, Block, Constraint, NonNegativeReals, summation, SolverFactory, minimize
def build_sp_model(yields):
'''
Rewritten version of the stochastic programming model using blocks.
Arguments:
yields: Yield information as a list, following the rank [wheat, corn, beets]
Return:
model: farmer problem model using blocks
'''
# Model
model = ConcreteModel()
# Sets
all_crops = ["WHEAT", "CORN", "BEETS"]
purchase_crops = ["WHEAT", "CORN"]
sell_crops = ["WHEAT", "CORN", "BEETS_FAVORABLE", "BEETS_UNFAVORABLE"]
scenarios = ["ABOVE", "AVERAGE", "BELOW"]
# First-stage decision: how much to land to dedicate to each crop
model.X = Var(all_crops, within=NonNegativeReals)
# Define a block for each scenario
def scenario_block_rule(b, scenario):
# Second-stage decision variables
b.Y = Var(purchase_crops, within=NonNegativeReals) # How much to purchase in this scenario
b.W = Var(sell_crops, within=NonNegativeReals) # How much to sell in this scenario
# Purchase cost and sales revenue for this scenario
if scenario == "ABOVE":
b.purchase_cost = 238 * b.Y["WHEAT"] + 210 * b.Y["CORN"]
b.sales_revenue = (
170 * b.W["WHEAT"] + 150 * b.W["CORN"]
+ 36 * b.W["BEETS_FAVORABLE"] + 10 * b.W["BEETS_UNFAVORABLE"]
)
elif scenario == "AVERAGE":
b.purchase_cost = 238 * b.Y["WHEAT"] + 210 * b.Y["CORN"]
b.sales_revenue = (
170 * b.W["WHEAT"] + 150 * b.W["CORN"]
+ 36 * b.W["BEETS_FAVORABLE"] + 10 * b.W["BEETS_UNFAVORABLE"]
)
else: # BELOW
b.purchase_cost = 238 * b.Y["WHEAT"] + 210 * b.Y["CORN"]
b.sales_revenue = (
170 * b.W["WHEAT"] + 150 * b.W["CORN"]
+ 36 * b.W["BEETS_FAVORABLE"] + 10 * b.W["BEETS_UNFAVORABLE"]
)
# Scenario constraints
if scenario == "ABOVE":
b.wheat_constraint = Constraint(expr=yields[0] * 1.2 * model.X["WHEAT"] + b.Y["WHEAT"] - b.W["WHEAT"] >= 200)
b.corn_constraint = Constraint(expr=yields[1] * 1.2 * model.X["CORN"] + b.Y["CORN"] - b.W["CORN"] >= 240)
b.beets_constraint = Constraint(expr=yields[2] * 1.2 * model.X["BEETS"]
- b.W["BEETS_FAVORABLE"] - b.W["BEETS_UNFAVORABLE"] >= 0)
elif scenario == "AVERAGE":
b.wheat_constraint = Constraint(expr=yields[0] * model.X["WHEAT"] + b.Y["WHEAT"] - b.W["WHEAT"] >= 200)
b.corn_constraint = Constraint(expr=yields[1] * model.X["CORN"] + b.Y["CORN"] - b.W["CORN"] >= 240)
b.beets_constraint = Constraint(expr=yields[2] * model.X["BEETS"]
- b.W["BEETS_FAVORABLE"] - b.W["BEETS_UNFAVORABLE"] >= 0)
else: # BELOW
b.wheat_constraint = Constraint(expr=yields[0] * 0.8 * model.X["WHEAT"] + b.Y["WHEAT"] - b.W["WHEAT"] >= 200)
b.corn_constraint = Constraint(expr=yields[1] * 0.8 * model.X["CORN"] + b.Y["CORN"] - b.W["CORN"] >= 240)
b.beets_constraint = Constraint(expr=yields[2] * 0.8 * model.X["BEETS"]
- b.W["BEETS_FAVORABLE"] - b.W["BEETS_UNFAVORABLE"] >= 0)
# Set upper bounds for BEETS_FAVORABLE
b.W["BEETS_FAVORABLE"].setub(6000)
# Create blocks for each scenario
model.scenarios = Block(scenarios, rule=scenario_block_rule)
# Objective function
def objective_rule(m):
planting_cost = 150 * m.X["WHEAT"] + 230 * m.X["CORN"] + 260 * m.X["BEETS"]
expected_purchase_cost = (1/3) * sum(m.scenarios[sc].purchase_cost for sc in scenarios)
expected_sales_revenue = (1/3) * sum(m.scenarios[sc].sales_revenue for sc in scenarios)
return planting_cost + expected_purchase_cost - expected_sales_revenue
model.OBJ = Objective(rule=objective_rule, sense=minimize)
# First-stage constraint: total area allocated to crops should not exceed 500
model.total_land_constraint = Constraint(expr=summation(model.X) <= 500)
return model
yields = [2.5, 3.0, 20.0]
model = build_sp_model(yields)
solver = SolverFactory('cbc')
solver.solve(model)
# Display the results
print("Planting decisions:")
for crop in model.X:
print(f"{crop}: {model.X[crop]()}")
Planting decisions:
WHEAT: 170.0
CORN: 80.0
BEETS: 250.0
We got the same answer as Stochastic Programming!
4.2.3. Seperate the Data from the Model#
The above blocks example is not really a big improvement from our alternate implementation in Stochastic Programming.
In the above code (from ChatGPT), we see:
Data for each scenario are hardcoded into the model
If statements are used to toggle the correct constraints for each scenario
The above code is not general purpose; updating it to use alternate scenario data or consider additional crops would take a lot of manual effort. It would be easy to make mistakes.
Thus, it is best practice to always seperate your specific problem data from the general mathematical model. Let’s see an example.
4.2.3.1. Pandas DataFrames#
Let’s use a pandas dataframe to store our problem data.
import pandas as pd
nominal_data = pd.read_csv("https://raw.githubusercontent.com/ndcbe/optimization/main/notebooks/data/farmers.csv")
nominal_data.head()
| Unnamed: 0 | Wheat | Corn | Beats | Units | |
|---|---|---|---|---|---|
| 0 | Yield | 2.5 | 3.0 | 20.0 | T/acre |
| 1 | Planting Cost | 150.0 | 230.0 | 260.0 | USD/acre |
| 2 | Favorable Selling Price | 170.0 | 150.0 | 36.0 | USD/T |
| 3 | Unfavorable Selling Price | NaN | NaN | 10.0 | USD/T |
| 4 | Purchase Price | 238.0 | 210.0 | NaN | USD/T |
This looks great. But it is best practice to have the rows be instances of data and the columns to be the types of data. We need to transpose the CSV file. Also, let’s drop the units for simplicity. ChatGPT is actually really helpful for transposing the CSV file!
nominal_data = pd.read_csv("https://raw.githubusercontent.com/ndcbe/optimization/main/notebooks/data/farmers2.csv")
nominal_data.head()
| index | Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | |
|---|---|---|---|---|---|---|---|---|
| 0 | Wheat | 2.5 | 150 | 170 | NaN | 238 | 200 | NaN |
| 1 | Corn | 3 | 230 | 150 | NaN | 210 | 240 | NaN |
| 2 | Beats | 20 | 260 | 36 | 10.0 | NaN | NaN | 6000 |
| 3 | Units | T/acre | USD/acre | USD/T | NaN | USD/T | T | T |
Now let’s drop the units for simplicity.
nominal_data = nominal_data.set_index("index")
nominal_data.head()
| Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | |
|---|---|---|---|---|---|---|---|
| index | |||||||
| Wheat | 2.5 | 150 | 170 | NaN | 238 | 200 | NaN |
| Corn | 3 | 230 | 150 | NaN | 210 | 240 | NaN |
| Beats | 20 | 260 | 36 | 10.0 | NaN | NaN | 6000 |
| Units | T/acre | USD/acre | USD/T | NaN | USD/T | T | T |
nominal_data.drop("Units", inplace=True)
nominal_data.head()
| Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | |
|---|---|---|---|---|---|---|---|
| index | |||||||
| Wheat | 2.5 | 150 | 170 | NaN | 238 | 200 | NaN |
| Corn | 3 | 230 | 150 | NaN | 210 | 240 | NaN |
| Beats | 20 | 260 | 36 | 10.0 | NaN | NaN | 6000 |
Finally, we need to convert the data entries from strings to numbers.
# Convert all elements to numbers
nominal_data = nominal_data.map(lambda x: pd.to_numeric(x, errors='coerce'))
nominal_data.head()
| Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | |
|---|---|---|---|---|---|---|---|
| index | |||||||
| Wheat | 2.5 | 150 | 170 | NaN | 238.0 | 200.0 | NaN |
| Corn | 3.0 | 230 | 150 | NaN | 210.0 | 240.0 | NaN |
| Beats | 20.0 | 260 | 36 | 10.0 | NaN | NaN | 6000.0 |
Finally, we should remove the NaN values and add columns with booleans.
nominal_data['Enforce Max Production'] = False
nominal_data['Enforce Min Requirement'] = False
nominal_data['Allow Purchases'] = True
nominal_data.head()
| Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | Enforce Max Production | Enforce Min Requirement | Allow Purchases | |
|---|---|---|---|---|---|---|---|---|---|---|
| index | ||||||||||
| Wheat | 2.5 | 150 | 170 | NaN | 238.0 | 200.0 | NaN | False | False | True |
| Corn | 3.0 | 230 | 150 | NaN | 210.0 | 240.0 | NaN | False | False | True |
| Beats | 20.0 | 260 | 36 | 10.0 | NaN | NaN | 6000.0 | False | False | True |
import numpy as np
for i, (index, row) in enumerate(nominal_data.iterrows()):
# print(row)
if not np.isnan(row['Maximum Favorable Production']):
nominal_data.at[index, 'Enforce Max Production'] = True
else:
nominal_data.at[index, 'Maximum Favorable Production'] = np.inf
nominal_data.at[index, 'Enforce Max Production'] = False
if not np.isnan(row['Minimum Requirement']):
nominal_data.at[index, 'Enforce Min Requirement'] = True
else:
nominal_data.at[index, 'Minimum Requirement'] = 0
nominal_data.at[index, 'Enforce Min Requirement'] = False
if np.isnan(row['Purchase Price']):
nominal_data.at[index, 'Allow Purchases'] = False
nominal_data.at[index, 'Purchase Price'] = 0
else:
nominal_data.at[index, 'Allow Purchases'] = True
if np.isnan(row['Unfavorable Selling Price']):
nominal_data.at[index, 'Unfavorable Selling Price'] = 0
nominal_data.head()
| Yield | Planting Cost | Favorable Selling Price | Unfavorable Selling Price | Purchase Price | Minimum Requirement | Maximum Favorable Production | Enforce Max Production | Enforce Min Requirement | Allow Purchases | |
|---|---|---|---|---|---|---|---|---|---|---|
| index | ||||||||||
| Wheat | 2.5 | 150 | 170 | 0.0 | 238.0 | 200.0 | inf | False | True | True |
| Corn | 3.0 | 230 | 150 | 0.0 | 210.0 | 240.0 | inf | False | True | True |
| Beats | 20.0 | 260 | 36 | 10.0 | 0.0 | 0.0 | 6000.0 | True | False | False |
4.2.3.2. Single Scenario Optimization Problem#
from pyomo.environ import Set, Param, Constraint, maximize
import numpy as np
def build_deterministic_model(nominal_data, m=None, skip_objective=False):
# Create concrete model (if not provided)
if m is None:
m = ConcreteModel()
# Sets
crops = nominal_data.index.to_list()
m.CROPS = Set(initialize=crops)
m.max_area = 500
# Parameters
m.cost = Param(m.CROPS, initialize=nominal_data["Planting Cost"], within=NonNegativeReals)
m.sell_price_favorable = Param(m.CROPS, initialize=nominal_data["Favorable Selling Price"], within=NonNegativeReals)
m.sell_price_unfavorable = Param(m.CROPS, initialize=nominal_data["Unfavorable Selling Price"], within=NonNegativeReals)
m.purchase_price = Param(m.CROPS, initialize=nominal_data["Purchase Price"], within=NonNegativeReals)
m.crop_yield = Param(m.CROPS, initialize=nominal_data["Yield"], within=NonNegativeReals)
m.min_required = Param(m.CROPS, initialize=nominal_data["Minimum Requirement"], within=NonNegativeReals)
m.max_possible = Param(m.CROPS, initialize=nominal_data["Maximum Favorable Production"], within=NonNegativeReals)
# Extract boolean parameters
m.enforce_max_production = nominal_data["Enforce Max Production"].to_dict()
m.enforce_min_requirement = nominal_data["Enforce Min Requirement"].to_dict()
m.allow_purchases = nominal_data["Allow Purchases"].to_dict()
# Stage 1 decision variables
m.plant = Var(m.CROPS, within=NonNegativeReals)
# Stage 2 decision variables
m.buy = Var(m.CROPS, within=NonNegativeReals) # purchases
m.sell_favorable = Var(m.CROPS, within=NonNegativeReals) # sales
m.sell_unfavorable = Var(m.CROPS, within=NonNegativeReals) # sales
for c in crops:
# Disable purchases for crops with NaN prices
if not m.allow_purchases[c]:
m.buy[c].fix(0)
# Enforce maximum production limits
if m.enforce_max_production[c]:
m.sell_favorable[c].setub(m.max_possible[c])
# Total area constraint
@m.Constraint()
def total_area(m):
return sum(m.plant[crop] for crop in m.CROPS) <= m.max_area
# Constraint on the production of each crop
@m.Constraint(m.CROPS)
def crop_min_constraint(m, crop):
if m.enforce_min_requirement[crop]:
return m.plant[crop] * m.crop_yield[crop] + m.buy[crop] - m.sell_unfavorable[crop] - m.sell_favorable[crop] >= m.min_required[crop]
else:
return m.plant[crop] * m.crop_yield[crop] + m.buy[crop] - m.sell_unfavorable[crop] - m.sell_favorable[crop] >= 0
# Maximize net profit
@m.Expression()
def net_profit(m):
return -sum(m.plant[c] * m.cost[c] for c in crops) - sum(m.buy[c] * m.purchase_price[c] for c in crops) + sum(m.sell_favorable[c] * m.sell_price_favorable[c] + m.sell_unfavorable[c] * m.sell_price_unfavorable[c] for c in crops)
if not skip_objective:
m.obj = Objective(expr=m.net_profit, sense=maximize)
return m
m = build_deterministic_model(nominal_data)
m.pprint()
1 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
7 Param Declarations
cost : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 20.0
Corn : 3.0
Wheat : 2.5
max_possible : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
4 Var Declarations
buy : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
plant : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_favorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : 6000.0 : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_unfavorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
1 Expression Declarations
net_profit : Size=1, Index=None
Key : Expression
None : - (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*buy[Wheat] + 210.0*buy[Corn] + 0.0*buy[Beats]) + (170*sell_favorable[Wheat] + 0.0*sell_unfavorable[Wheat] + 150*sell_favorable[Corn] + 0.0*sell_unfavorable[Corn] + 36*sell_favorable[Beats] + 10.0*sell_unfavorable[Beats])
1 Objective Declarations
obj : Size=1, Index=None, Active=True
Key : Active : Sense : Expression
None : True : maximize : net_profit
2 Constraint Declarations
crop_min_constraint : Size=3, Index=CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 20.0*plant[Beats] + buy[Beats] - sell_unfavorable[Beats] - sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 3.0*plant[Corn] + buy[Corn] - sell_unfavorable[Corn] - sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 2.5*plant[Wheat] + buy[Wheat] - sell_unfavorable[Wheat] - sell_favorable[Wheat] : +Inf : True
total_area : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : plant[Wheat] + plant[Corn] + plant[Beats] : 500.0 : True
16 Declarations: CROPS cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible plant buy sell_favorable sell_unfavorable total_area crop_min_constraint net_profit obj
solver = SolverFactory('ipopt')
solver.solve(m, tee=True)
Ipopt 3.13.2:
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit http://projects.coin-or.org/Ipopt
This version of Ipopt was compiled from source code available at
https://github.com/IDAES/Ipopt as part of the Institute for the Design of
Advanced Energy Systems Process Systems Engineering Framework (IDAES PSE
Framework) Copyright (c) 2018-2019. See https://github.com/IDAES/idaes-pse.
This version of Ipopt was compiled using HSL, a collection of Fortran codes
for large-scale scientific computation. All technical papers, sales and
publicity material resulting from use of the HSL codes within IPOPT must
contain the following acknowledgement:
HSL, a collection of Fortran codes for large-scale scientific
computation. See http://www.hsl.rl.ac.uk.
******************************************************************************
This is Ipopt version 3.13.2, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 0
Number of nonzeros in inequality constraint Jacobian.: 14
Number of nonzeros in Lagrangian Hessian.............: 0
Total number of variables............................: 11
variables with only lower bounds: 10
variables with lower and upper bounds: 1
variables with only upper bounds: 0
Total number of equality constraints.................: 0
Total number of inequality constraints...............: 4
inequality constraints with only lower bounds: 3
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 7.2199928e+00 2.40e+02 5.79e+01 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 1.6884409e+01 2.40e+02 5.72e+01 -1.0 2.34e+02 - 4.81e-03 1.01e-02h 1
2 3.9285122e+02 2.37e+02 5.62e+01 -1.0 1.59e+02 - 3.80e-04 1.24e-02h 1
3 8.3995684e+02 2.33e+02 5.45e+01 -1.0 1.59e+02 - 2.68e-02 1.51e-02h 1
4 3.2827735e+03 2.11e+02 5.56e+01 -1.0 1.70e+02 - 2.75e-04 9.58e-02h 1
5 -1.2649345e+05 2.09e+02 5.52e+01 -1.0 1.10e+05 - 3.11e-05 1.07e-02f 1
6 -1.2720560e+05 2.07e+02 5.18e+01 -1.0 9.48e+02 - 1.36e-01 9.17e-03f 1
7 -1.1792715e+05 1.77e+02 4.66e+01 -1.0 3.44e+02 - 3.17e-02 1.45e-01h 1
8 -9.9006212e+04 1.15e+02 3.69e+01 -1.0 3.09e+02 - 8.94e-03 3.52e-01h 1
9 -1.4404818e+05 8.39e+01 4.94e+01 -1.0 2.20e+04 - 1.93e-03 2.68e-01f 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 -1.3875045e+05 6.50e+01 2.90e+01 -1.0 2.62e+02 - 1.00e+00 2.25e-01h 1
11 -1.1859714e+05 0.00e+00 1.00e-06 -1.0 1.08e+02 - 1.00e+00 1.00e+00h 1
12 -1.1859992e+05 0.00e+00 2.83e-08 -2.5 2.73e-02 - 1.00e+00 1.00e+00f 1
13 -1.1860000e+05 0.00e+00 1.50e-09 -3.8 7.52e-04 - 1.00e+00 1.00e+00f 1
14 -1.1860000e+05 0.00e+00 1.85e-11 -5.7 4.18e-05 - 1.00e+00 1.00e+00f 1
15 -1.1860000e+05 0.00e+00 3.52e-14 -8.6 5.18e-07 - 1.00e+00 1.00e+00f 1
Number of Iterations....: 15
(scaled) (unscaled)
Objective...............: -4.5615385645779483e+04 -1.1860000267902664e+05
Dual infeasibility......: 3.5194748420182662e-14 9.1506345892474909e-14
Constraint violation....: 0.0000000000000000e+00 0.0000000000000000e+00
Complementarity.........: 2.5071227145647319e-09 6.5185190578683025e-09
Overall NLP error.......: 2.5071227145647319e-09 6.5185190578683025e-09
Number of objective function evaluations = 16
Number of objective gradient evaluations = 16
Number of equality constraint evaluations = 0
Number of inequality constraint evaluations = 16
Number of equality constraint Jacobian evaluations = 0
Number of inequality constraint Jacobian evaluations = 16
Number of Lagrangian Hessian evaluations = 15
Total CPU secs in IPOPT (w/o function evaluations) = 0.002
Total CPU secs in NLP function evaluations = 0.000
EXIT: Optimal Solution Found.
{'Problem': [{'Lower bound': -inf, 'Upper bound': inf, 'Number of objectives': 1, 'Number of constraints': 4, 'Number of variables': 11, 'Sense': 'unknown'}], 'Solver': [{'Status': 'ok', 'Message': 'Ipopt 3.13.2\\x3a Optimal Solution Found', 'Termination condition': 'optimal', 'Id': 0, 'Error rc': 0, 'Time': 0.021295785903930664}], 'Solution': [OrderedDict([('number of solutions', 0), ('number of solutions displayed', 0)])]}
def print_solution(m):
# Display the results
print("Planting decisions:")
for crop in m.CROPS:
print(f"{crop}: {round(m.plant[crop](),2)} acres")
print("\nPurchase decisions:")
for crop in m.CROPS:
print(f"{crop}: {round(m.buy[crop](),2)} T")
print("\nSales decisions (favorable price):")
for crop in m.CROPS:
print(f"{crop}: {round(m.sell_favorable[crop](),2)} T")
print("\nSales decisions (unfavorable price):")
for crop in m.CROPS:
print(f"{crop}: {round(m.sell_unfavorable[crop](),2)} T")
print_solution(m)
Planting decisions:
Wheat: 120.0 acres
Corn: 80.0 acres
Beats: 300.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 100.0 T
Corn: 0.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
4.2.4. Blocks#
Now we will use blocks to construct the two-stage stochastic program.
4.2.4.1. Prepare Input Data#
scenarios = {}
yield_factors = [0.8, 1.0, 1.2]
scienario_names = ["BELOW", "AVERAGE", "ABOVE"]
for i, name in enumerate(scienario_names):
# Copy the data
scenario = nominal_data.copy()
# Change the yield
scenario["Yield"] = scenario["Yield"] * yield_factors[i]
scenarios[name] = scenario
print(scenarios)
{'BELOW': Yield Planting Cost Favorable Selling Price \
index
Wheat 2.0 150 170
Corn 2.4 230 150
Beats 16.0 260 36
Unfavorable Selling Price Purchase Price Minimum Requirement \
index
Wheat 0.0 238.0 200.0
Corn 0.0 210.0 240.0
Beats 10.0 0.0 0.0
Maximum Favorable Production Enforce Max Production \
index
Wheat inf False
Corn inf False
Beats 6000.0 True
Enforce Min Requirement Allow Purchases
index
Wheat True True
Corn True True
Beats False False , 'AVERAGE': Yield Planting Cost Favorable Selling Price \
index
Wheat 2.5 150 170
Corn 3.0 230 150
Beats 20.0 260 36
Unfavorable Selling Price Purchase Price Minimum Requirement \
index
Wheat 0.0 238.0 200.0
Corn 0.0 210.0 240.0
Beats 10.0 0.0 0.0
Maximum Favorable Production Enforce Max Production \
index
Wheat inf False
Corn inf False
Beats 6000.0 True
Enforce Min Requirement Allow Purchases
index
Wheat True True
Corn True True
Beats False False , 'ABOVE': Yield Planting Cost Favorable Selling Price \
index
Wheat 3.0 150 170
Corn 3.6 230 150
Beats 24.0 260 36
Unfavorable Selling Price Purchase Price Minimum Requirement \
index
Wheat 0.0 238.0 200.0
Corn 0.0 210.0 240.0
Beats 10.0 0.0 0.0
Maximum Favorable Production Enforce Max Production \
index
Wheat inf False
Corn inf False
Beats 6000.0 True
Enforce Min Requirement Allow Purchases
index
Wheat True True
Corn True True
Beats False False }
4.2.4.2. Solve the Deterministic Model for Each Stage#
for i, (key, value) in enumerate(scenarios.items()):
print(f"*** Solving scenario {i} = {key} ***")
m = build_deterministic_model(value)
solver = SolverFactory('ipopt')
solver.solve(m, tee=False)
print_solution(m)
print("\n")
*** Solving scenario 0 = BELOW ***
Planting decisions:
Wheat: 100.0 acres
Corn: 25.0 acres
Beats: 375.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 180.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Solving scenario 1 = AVERAGE ***
Planting decisions:
Wheat: 120.0 acres
Corn: 80.0 acres
Beats: 300.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 100.0 T
Corn: 0.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Solving scenario 2 = ABOVE ***
Planting decisions:
Wheat: 183.33 acres
Corn: 66.67 acres
Beats: 250.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 350.0 T
Corn: 0.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
4.2.4.3. Define the Two-Stage Stochastic Program Using Blocks#
from pyomo.environ import Expression
def create_two_stage_stochastic_model(scenarios):
# Model
m = ConcreteModel()
nominal_scenario_name = "AVERAGE"
nominal_data = scenarios[nominal_scenario_name]
# Sets
crops = nominal_data.index.to_list()
m.CROPS = Set(initialize=crops)
m.max_area = 500
# Stage 1 variables
m.plant = Var(m.CROPS, within=NonNegativeReals, initialize=100, bounds=(0, m.max_area))
# Stage 1 constraint
@m.Constraint()
def total_area(m):
return sum(m.plant[crop] for crop in m.CROPS) <= m.max_area
def second_stage_block(b, scenario_name):
# print("name = ", scenario_name)
# Grab data for this specific scenario
s_data = scenarios[scenario_name]
# Parameters
b.cost = Param(m.CROPS, initialize=s_data["Planting Cost"], within=NonNegativeReals)
b.sell_price_favorable = Param(m.CROPS, initialize=s_data["Favorable Selling Price"], within=NonNegativeReals)
b.sell_price_unfavorable = Param(m.CROPS, initialize=s_data["Unfavorable Selling Price"], within=NonNegativeReals)
b.purchase_price = Param(m.CROPS, initialize=s_data["Purchase Price"], within=NonNegativeReals)
b.crop_yield = Param(m.CROPS, initialize=s_data["Yield"], within=NonNegativeReals)
b.min_required = Param(m.CROPS, initialize=s_data["Minimum Requirement"], within=NonNegativeReals)
b.max_possible = Param(m.CROPS, initialize=s_data["Maximum Favorable Production"], within=NonNegativeReals)
# Extract boolean parameters
b.enforce_max_production = s_data["Enforce Max Production"].to_dict()
b.enforce_min_requirement = s_data["Enforce Min Requirement"].to_dict()
b.allow_purchases = s_data["Allow Purchases"].to_dict()
# Stage 2 decision variables
b.buy = Var(m.CROPS, within=NonNegativeReals, initialize=0) # purchases
b.sell_favorable = Var(m.CROPS, within=NonNegativeReals, initialize=0) # sales
b.sell_unfavorable = Var(m.CROPS, within=NonNegativeReals, initialize=0) # sales
for c in crops:
# Disable purchases for crops with NaN prices
if not b.allow_purchases[c]:
b.buy[c].fix(0)
# Enforce maximum production limits
if b.enforce_max_production[c]:
b.sell_favorable[c].setub(b.max_possible[c])
# Constraint on the production of each crop
@b.Constraint(m.CROPS)
def crop_min_constraint(b, crop):
if b.enforce_min_requirement[crop]:
return m.plant[crop] * b.crop_yield[crop] + b.buy[crop] - b.sell_unfavorable[crop] - b.sell_favorable[crop] >= b.min_required[crop]
else:
return m.plant[crop] * b.crop_yield[crop] + b.buy[crop] - b.sell_unfavorable[crop] - b.sell_favorable[crop] >= 0
@b.Expression()
def scenario_net_profit(b):
return -sum(m.plant[c] * b.cost[c] for c in crops) - sum(b.buy[c] * b.purchase_price[c] for c in crops) + sum(b.sell_favorable[c] * b.sell_price_favorable[c] + b.sell_unfavorable[c] * b.sell_price_unfavorable[c] for c in crops)
# Create blocks for each scenario
m.planning_scenarios = Block(scenarios.keys(), rule=second_stage_block)
# Maximize net profot
@m.Objective(sense=maximize)
def obj(m):
return sum(m.planning_scenarios[sc].scenario_net_profit for sc in scenarios.keys())/len(scenarios)
return m
m = create_two_stage_stochastic_model(scenarios)
m.pprint()
1 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
1 Var Declarations
plant : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 100 : 500 : False : True : NonNegativeReals
Corn : 0 : 100 : 500 : False : True : NonNegativeReals
Wheat : 0 : 100 : 500 : False : True : NonNegativeReals
1 Objective Declarations
obj : Size=1, Index=None, Active=True
Key : Active : Sense : Expression
None : True : maximize : ((- (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[BELOW].buy[Wheat] + 210.0*planning_scenarios[BELOW].buy[Corn] + 0.0*planning_scenarios[BELOW].buy[Beats]) + (170*planning_scenarios[BELOW].sell_favorable[Wheat] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Wheat] + 150*planning_scenarios[BELOW].sell_favorable[Corn] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Corn] + 36*planning_scenarios[BELOW].sell_favorable[Beats] + 10.0*planning_scenarios[BELOW].sell_unfavorable[Beats])) + (- (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[AVERAGE].buy[Wheat] + 210.0*planning_scenarios[AVERAGE].buy[Corn] + 0.0*planning_scenarios[AVERAGE].buy[Beats]) + (170*planning_scenarios[AVERAGE].sell_favorable[Wheat] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Wheat] + 150*planning_scenarios[AVERAGE].sell_favorable[Corn] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Corn] + 36*planning_scenarios[AVERAGE].sell_favorable[Beats] + 10.0*planning_scenarios[AVERAGE].sell_unfavorable[Beats])) + (- (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[ABOVE].buy[Wheat] + 210.0*planning_scenarios[ABOVE].buy[Corn] + 0.0*planning_scenarios[ABOVE].buy[Beats]) + (170*planning_scenarios[ABOVE].sell_favorable[Wheat] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Wheat] + 150*planning_scenarios[ABOVE].sell_favorable[Corn] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Corn] + 36*planning_scenarios[ABOVE].sell_favorable[Beats] + 10.0*planning_scenarios[ABOVE].sell_unfavorable[Beats])))/3
1 Constraint Declarations
total_area : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : plant[Wheat] + plant[Corn] + plant[Beats] : 500.0 : True
1 Block Declarations
planning_scenarios : Size=3, Index={BELOW, AVERAGE, ABOVE}, Active=True
planning_scenarios[ABOVE] : Active=True
7 Param Declarations
cost : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 24.0
Corn : 3.5999999999999996
Wheat : 3.0
max_possible : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
3 Var Declarations
buy : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
sell_favorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : 6000.0 : False : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
sell_unfavorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : False : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
1 Expression Declarations
scenario_net_profit : Size=1, Index=None
Key : Expression
None : - (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[ABOVE].buy[Wheat] + 210.0*planning_scenarios[ABOVE].buy[Corn] + 0.0*planning_scenarios[ABOVE].buy[Beats]) + (170*planning_scenarios[ABOVE].sell_favorable[Wheat] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Wheat] + 150*planning_scenarios[ABOVE].sell_favorable[Corn] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Corn] + 36*planning_scenarios[ABOVE].sell_favorable[Beats] + 10.0*planning_scenarios[ABOVE].sell_unfavorable[Beats])
1 Constraint Declarations
crop_min_constraint : Size=3, Index=CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 24.0*plant[Beats] + planning_scenarios[ABOVE].buy[Beats] - planning_scenarios[ABOVE].sell_unfavorable[Beats] - planning_scenarios[ABOVE].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 3.5999999999999996*plant[Corn] + planning_scenarios[ABOVE].buy[Corn] - planning_scenarios[ABOVE].sell_unfavorable[Corn] - planning_scenarios[ABOVE].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 3.0*plant[Wheat] + planning_scenarios[ABOVE].buy[Wheat] - planning_scenarios[ABOVE].sell_unfavorable[Wheat] - planning_scenarios[ABOVE].sell_favorable[Wheat] : +Inf : True
12 Declarations: cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible buy sell_favorable sell_unfavorable crop_min_constraint scenario_net_profit
planning_scenarios[AVERAGE] : Active=True
7 Param Declarations
cost : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 20.0
Corn : 3.0
Wheat : 2.5
max_possible : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
3 Var Declarations
buy : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
sell_favorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : 6000.0 : False : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
sell_unfavorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : False : False : NonNegativeReals
Corn : 0 : 0 : None : False : False : NonNegativeReals
Wheat : 0 : 0 : None : False : False : NonNegativeReals
1 Expression Declarations
scenario_net_profit : Size=1, Index=None
Key : Expression
None : - (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[AVERAGE].buy[Wheat] + 210.0*planning_scenarios[AVERAGE].buy[Corn] + 0.0*planning_scenarios[AVERAGE].buy[Beats]) + (170*planning_scenarios[AVERAGE].sell_favorable[Wheat] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Wheat] + 150*planning_scenarios[AVERAGE].sell_favorable[Corn] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Corn] + 36*planning_scenarios[AVERAGE].sell_favorable[Beats] + 10.0*planning_scenarios[AVERAGE].sell_unfavorable[Beats])
1 Constraint Declarations
crop_min_constraint : Size=3, Index=CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 20.0*plant[Beats] + planning_scenarios[AVERAGE].buy[Beats] - planning_scenarios[AVERAGE].sell_unfavorable[Beats] - planning_scenarios[AVERAGE].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 3.0*plant[Corn] + planning_scenarios[AVERAGE].buy[Corn] - planning_scenarios[AVERAGE].sell_unfavorable[Corn] - planning_scenarios[AVERAGE].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 2.5*plant[Wheat] + planning_scenarios[AVERAGE].buy[Wheat] - planning_scenarios[AVERAGE].sell_unfavorable[Wheat] - planning_scenarios[AVERAGE].sell_favorable[Wheat] : +Inf : True
12 Declarations: cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible buy sell_favorable sell_unfavorable crop_min_constraint scenario_net_profit
planning_scenarios[BELOW] : Active=True
7 Param Declarations
cost : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 16.0
Corn : 2.4000000000000004
Wheat : 2.0
max_possible : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
3 Var Declarations
buy : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : 0 : None : False : True : NonNegativeReals
Wheat : 0 : 0 : None : False : True : NonNegativeReals
sell_favorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : 6000.0 : False : True : NonNegativeReals
Corn : 0 : 0 : None : False : True : NonNegativeReals
Wheat : 0 : 0 : None : False : True : NonNegativeReals
sell_unfavorable : Size=3, Index=CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : False : True : NonNegativeReals
Corn : 0 : 0 : None : False : True : NonNegativeReals
Wheat : 0 : 0 : None : False : True : NonNegativeReals
1 Expression Declarations
scenario_net_profit : Size=1, Index=None
Key : Expression
None : - (150*plant[Wheat] + 230*plant[Corn] + 260*plant[Beats]) - (238.0*planning_scenarios[BELOW].buy[Wheat] + 210.0*planning_scenarios[BELOW].buy[Corn] + 0.0*planning_scenarios[BELOW].buy[Beats]) + (170*planning_scenarios[BELOW].sell_favorable[Wheat] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Wheat] + 150*planning_scenarios[BELOW].sell_favorable[Corn] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Corn] + 36*planning_scenarios[BELOW].sell_favorable[Beats] + 10.0*planning_scenarios[BELOW].sell_unfavorable[Beats])
1 Constraint Declarations
crop_min_constraint : Size=3, Index=CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 16.0*plant[Beats] + planning_scenarios[BELOW].buy[Beats] - planning_scenarios[BELOW].sell_unfavorable[Beats] - planning_scenarios[BELOW].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 2.4000000000000004*plant[Corn] + planning_scenarios[BELOW].buy[Corn] - planning_scenarios[BELOW].sell_unfavorable[Corn] - planning_scenarios[BELOW].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 2.0*plant[Wheat] + planning_scenarios[BELOW].buy[Wheat] - planning_scenarios[BELOW].sell_unfavorable[Wheat] - planning_scenarios[BELOW].sell_favorable[Wheat] : +Inf : True
12 Declarations: cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible buy sell_favorable sell_unfavorable crop_min_constraint scenario_net_profit
5 Declarations: CROPS plant total_area planning_scenarios obj
4.2.4.4. Solve and Inspect Solution#
solver = SolverFactory('ipopt')
solver.solve(m, tee=True)
Ipopt 3.13.2:
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit http://projects.coin-or.org/Ipopt
This version of Ipopt was compiled from source code available at
https://github.com/IDAES/Ipopt as part of the Institute for the Design of
Advanced Energy Systems Process Systems Engineering Framework (IDAES PSE
Framework) Copyright (c) 2018-2019. See https://github.com/IDAES/idaes-pse.
This version of Ipopt was compiled using HSL, a collection of Fortran codes
for large-scale scientific computation. All technical papers, sales and
publicity material resulting from use of the HSL codes within IPOPT must
contain the following acknowledgement:
HSL, a collection of Fortran codes for large-scale scientific
computation. See http://www.hsl.rl.ac.uk.
******************************************************************************
This is Ipopt version 3.13.2, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 0
Number of nonzeros in inequality constraint Jacobian.: 36
Number of nonzeros in Lagrangian Hessian.............: 0
Total number of variables............................: 27
variables with only lower bounds: 21
variables with lower and upper bounds: 6
variables with only upper bounds: 0
Total number of equality constraints.................: 0
Total number of inequality constraints...............: 10
inequality constraints with only lower bounds: 9
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 6.4000820e+04 1.00e-02 2.05e+01 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 6.2395224e+04 2.31e+00 1.95e+01 -1.0 3.88e+03 - 4.49e-02 3.02e-02f 1
2 5.9184223e+04 2.16e+00 2.01e+01 -1.0 3.60e+03 - 4.66e-04 7.42e-02f 1
3 5.4100646e+04 2.10e+00 1.97e+01 -1.0 3.19e+03 - 2.19e-03 2.74e-02f 1
4 4.7681497e+04 2.00e+00 1.96e+01 -1.0 3.05e+03 - 4.94e-03 4.70e-02f 1
5 4.4800643e+04 1.97e+00 1.98e+01 -1.0 2.82e+03 - 7.76e-03 1.70e-02f 1
6 3.7505505e+04 1.85e+00 2.18e+01 -1.0 2.73e+03 - 6.44e-03 5.89e-02f 1
7 -2.9487452e+04 1.76e+00 2.17e+01 -1.0 1.51e+04 - 2.42e-04 4.73e-02f 1
8 -3.0999686e+04 1.74e+00 2.89e+01 -1.0 4.56e+03 - 3.56e-01 8.57e-03f 1
9 -5.1611237e+04 1.11e+00 1.38e+01 -1.0 2.54e+03 - 5.74e-03 3.63e-01f 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 -8.8197316e+04 5.99e-01 3.62e+01 -1.0 7.42e+03 - 1.72e-03 4.60e-01f 1
11 -1.0802365e+05 3.55e-01 2.64e+01 -1.0 6.28e+03 - 3.36e-01 4.06e-01f 1
12 -1.0827120e+05 2.46e-01 3.16e+01 -1.0 8.43e+01 - 8.82e-03 3.02e-01f 1
13 -1.0837183e+05 2.35e-01 1.82e+01 -1.0 1.18e+03 - 3.19e-01 4.36e-02f 1
14 -1.0837336e+05 2.26e-01 1.13e+01 -1.0 1.31e+01 - 9.83e-01 3.91e-02f 1
15 -1.0838298e+05 0.00e+00 2.19e-03 -1.0 5.22e+00 - 9.96e-01 1.00e+00f 1
16 -1.0838980e+05 0.00e+00 2.83e-08 -2.5 2.57e-01 - 1.00e+00 1.00e+00f 1
17 -1.0838999e+05 0.00e+00 1.50e-09 -3.8 7.44e-03 - 1.00e+00 1.00e+00f 1
18 -1.0839000e+05 0.00e+00 1.85e-11 -5.7 4.03e-04 - 1.00e+00 1.00e+00f 1
19 -1.0839000e+05 0.00e+00 2.67e-14 -8.6 5.00e-06 - 1.00e+00 1.00e+00f 1
Number of Iterations....: 19
(scaled) (unscaled)
Objective...............: -4.1688462536958956e+04 -1.0839000259609328e+05
Dual infeasibility......: 2.6663000555547025e-14 6.9323801444422252e-14
Constraint violation....: 0.0000000000000000e+00 0.0000000000000000e+00
Complementarity.........: 2.5071227145650024e-09 6.5185190578690056e-09
Overall NLP error.......: 2.5071227145650024e-09 6.5185190578690056e-09
Number of objective function evaluations = 20
Number of objective gradient evaluations = 20
Number of equality constraint evaluations = 0
Number of inequality constraint evaluations = 20
Number of equality constraint Jacobian evaluations = 0
Number of inequality constraint Jacobian evaluations = 20
Number of Lagrangian Hessian evaluations = 19
Total CPU secs in IPOPT (w/o function evaluations) = 0.002
Total CPU secs in NLP function evaluations = 0.000
EXIT: Optimal Solution Found.
{'Problem': [{'Lower bound': -inf, 'Upper bound': inf, 'Number of objectives': 1, 'Number of constraints': 10, 'Number of variables': 27, 'Sense': 'unknown'}], 'Solver': [{'Status': 'ok', 'Message': 'Ipopt 3.13.2\\x3a Optimal Solution Found', 'Termination condition': 'optimal', 'Id': 0, 'Error rc': 0, 'Time': 0.02140498161315918}], 'Solution': [OrderedDict([('number of solutions', 0), ('number of solutions displayed', 0)])]}
# Display the results
print("Planting decisions:")
for crop in m.CROPS:
print(f"{crop}: {round(m.plant[crop](),2)} acres")
for i, name in enumerate(scenarios.keys()):
print(f"\n *** Scenario {i} = {name} ***")
print("\n\tPurchase decisions:")
for crop in m.CROPS:
print(f"\t{crop}: {round(m.planning_scenarios[name].buy[crop](),2)} T")
print("\n\tSales decisions (favorable price):")
for crop in m.CROPS:
print(f"\t{crop}: {round(m.planning_scenarios[name].sell_favorable[crop](),2)} T")
print("\n\tSales decisions (unfavorable price):")
for crop in m.CROPS:
print(f"\t{crop}: {round(m.planning_scenarios[name].sell_unfavorable[crop](),2)} T")
Planting decisions:
Wheat: 170.0 acres
Corn: 80.0 acres
Beats: 250.0 acres
*** Scenario 0 = BELOW ***
Purchase decisions:
Wheat: 0.0 T
Corn: 48.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 140.0 T
Corn: 0.0 T
Beats: 4000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Scenario 1 = AVERAGE ***
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 225.0 T
Corn: 0.0 T
Beats: 5000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Scenario 2 = ABOVE ***
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 310.0 T
Corn: 48.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
4.2.5. Blocks with More Code Reuse#
How can we compactly write our problem reusing our function for the deterministic (single scenario) case?
def create_two_stage_stochastic_model_better(scenarios):
# Create a concrete model
m = ConcreteModel()
# Create a block for each scenario
# This creates a copy of your deterministic model for each scenario
m.planning_scenarios = Block(scenarios.keys(), rule=lambda b, s: build_deterministic_model(scenarios[s], m=b, skip_objective=True))
# Enforce the same planting (stage 1) decisions for all scenarios
m.SCENARIOS = Set(initialize=scenarios.keys())
m.CROPS = Set(initialize=scenarios["AVERAGE"].index.to_list())
@m.Constraint(m.SCENARIOS, m.CROPS)
def enforce_same_planting(m, s, c):
if s == "AVERAGE":
return Constraint.Skip
else:
return m.planning_scenarios[s].plant[c] == m.planning_scenarios["AVERAGE"].plant[c]
@m.Objective(sense=maximize)
def obj(m):
return sum(m.planning_scenarios[s].net_profit for s in scenarios.keys())/len(scenarios)
return m
m = create_two_stage_stochastic_model_better(scenarios)
m.pprint()
2 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
SCENARIOS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'BELOW', 'AVERAGE', 'ABOVE'}
1 Objective Declarations
obj : Size=1, Index=None, Active=True
Key : Active : Sense : Expression
None : True : maximize : ((- (150*planning_scenarios[BELOW].plant[Wheat] + 230*planning_scenarios[BELOW].plant[Corn] + 260*planning_scenarios[BELOW].plant[Beats]) - (238.0*planning_scenarios[BELOW].buy[Wheat] + 210.0*planning_scenarios[BELOW].buy[Corn] + 0.0*planning_scenarios[BELOW].buy[Beats]) + (170*planning_scenarios[BELOW].sell_favorable[Wheat] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Wheat] + 150*planning_scenarios[BELOW].sell_favorable[Corn] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Corn] + 36*planning_scenarios[BELOW].sell_favorable[Beats] + 10.0*planning_scenarios[BELOW].sell_unfavorable[Beats])) + (- (150*planning_scenarios[AVERAGE].plant[Wheat] + 230*planning_scenarios[AVERAGE].plant[Corn] + 260*planning_scenarios[AVERAGE].plant[Beats]) - (238.0*planning_scenarios[AVERAGE].buy[Wheat] + 210.0*planning_scenarios[AVERAGE].buy[Corn] + 0.0*planning_scenarios[AVERAGE].buy[Beats]) + (170*planning_scenarios[AVERAGE].sell_favorable[Wheat] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Wheat] + 150*planning_scenarios[AVERAGE].sell_favorable[Corn] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Corn] + 36*planning_scenarios[AVERAGE].sell_favorable[Beats] + 10.0*planning_scenarios[AVERAGE].sell_unfavorable[Beats])) + (- (150*planning_scenarios[ABOVE].plant[Wheat] + 230*planning_scenarios[ABOVE].plant[Corn] + 260*planning_scenarios[ABOVE].plant[Beats]) - (238.0*planning_scenarios[ABOVE].buy[Wheat] + 210.0*planning_scenarios[ABOVE].buy[Corn] + 0.0*planning_scenarios[ABOVE].buy[Beats]) + (170*planning_scenarios[ABOVE].sell_favorable[Wheat] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Wheat] + 150*planning_scenarios[ABOVE].sell_favorable[Corn] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Corn] + 36*planning_scenarios[ABOVE].sell_favorable[Beats] + 10.0*planning_scenarios[ABOVE].sell_unfavorable[Beats])))/3
1 Constraint Declarations
enforce_same_planting : Size=6, Index=SCENARIOS*CROPS, Active=True
Key : Lower : Body : Upper : Active
('ABOVE', 'Beats') : 0.0 : planning_scenarios[ABOVE].plant[Beats] - planning_scenarios[AVERAGE].plant[Beats] : 0.0 : True
('ABOVE', 'Corn') : 0.0 : planning_scenarios[ABOVE].plant[Corn] - planning_scenarios[AVERAGE].plant[Corn] : 0.0 : True
('ABOVE', 'Wheat') : 0.0 : planning_scenarios[ABOVE].plant[Wheat] - planning_scenarios[AVERAGE].plant[Wheat] : 0.0 : True
('BELOW', 'Beats') : 0.0 : planning_scenarios[BELOW].plant[Beats] - planning_scenarios[AVERAGE].plant[Beats] : 0.0 : True
('BELOW', 'Corn') : 0.0 : planning_scenarios[BELOW].plant[Corn] - planning_scenarios[AVERAGE].plant[Corn] : 0.0 : True
('BELOW', 'Wheat') : 0.0 : planning_scenarios[BELOW].plant[Wheat] - planning_scenarios[AVERAGE].plant[Wheat] : 0.0 : True
1 Block Declarations
planning_scenarios : Size=3, Index={BELOW, AVERAGE, ABOVE}, Active=True
planning_scenarios[ABOVE] : Active=True
1 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
7 Param Declarations
cost : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 24.0
Corn : 3.5999999999999996
Wheat : 3.0
max_possible : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=planning_scenarios[ABOVE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
4 Var Declarations
buy : Size=3, Index=planning_scenarios[ABOVE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
plant : Size=3, Index=planning_scenarios[ABOVE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_favorable : Size=3, Index=planning_scenarios[ABOVE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : 6000.0 : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_unfavorable : Size=3, Index=planning_scenarios[ABOVE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
1 Expression Declarations
net_profit : Size=1, Index=None
Key : Expression
None : - (150*planning_scenarios[ABOVE].plant[Wheat] + 230*planning_scenarios[ABOVE].plant[Corn] + 260*planning_scenarios[ABOVE].plant[Beats]) - (238.0*planning_scenarios[ABOVE].buy[Wheat] + 210.0*planning_scenarios[ABOVE].buy[Corn] + 0.0*planning_scenarios[ABOVE].buy[Beats]) + (170*planning_scenarios[ABOVE].sell_favorable[Wheat] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Wheat] + 150*planning_scenarios[ABOVE].sell_favorable[Corn] + 0.0*planning_scenarios[ABOVE].sell_unfavorable[Corn] + 36*planning_scenarios[ABOVE].sell_favorable[Beats] + 10.0*planning_scenarios[ABOVE].sell_unfavorable[Beats])
2 Constraint Declarations
crop_min_constraint : Size=3, Index=planning_scenarios[ABOVE].CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 24.0*planning_scenarios[ABOVE].plant[Beats] + planning_scenarios[ABOVE].buy[Beats] - planning_scenarios[ABOVE].sell_unfavorable[Beats] - planning_scenarios[ABOVE].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 3.5999999999999996*planning_scenarios[ABOVE].plant[Corn] + planning_scenarios[ABOVE].buy[Corn] - planning_scenarios[ABOVE].sell_unfavorable[Corn] - planning_scenarios[ABOVE].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 3.0*planning_scenarios[ABOVE].plant[Wheat] + planning_scenarios[ABOVE].buy[Wheat] - planning_scenarios[ABOVE].sell_unfavorable[Wheat] - planning_scenarios[ABOVE].sell_favorable[Wheat] : +Inf : True
total_area : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : planning_scenarios[ABOVE].plant[Wheat] + planning_scenarios[ABOVE].plant[Corn] + planning_scenarios[ABOVE].plant[Beats] : 500.0 : True
15 Declarations: CROPS cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible plant buy sell_favorable sell_unfavorable total_area crop_min_constraint net_profit
planning_scenarios[AVERAGE] : Active=True
1 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
7 Param Declarations
cost : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 20.0
Corn : 3.0
Wheat : 2.5
max_possible : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
4 Var Declarations
buy : Size=3, Index=planning_scenarios[AVERAGE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
plant : Size=3, Index=planning_scenarios[AVERAGE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_favorable : Size=3, Index=planning_scenarios[AVERAGE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : 6000.0 : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_unfavorable : Size=3, Index=planning_scenarios[AVERAGE].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
1 Expression Declarations
net_profit : Size=1, Index=None
Key : Expression
None : - (150*planning_scenarios[AVERAGE].plant[Wheat] + 230*planning_scenarios[AVERAGE].plant[Corn] + 260*planning_scenarios[AVERAGE].plant[Beats]) - (238.0*planning_scenarios[AVERAGE].buy[Wheat] + 210.0*planning_scenarios[AVERAGE].buy[Corn] + 0.0*planning_scenarios[AVERAGE].buy[Beats]) + (170*planning_scenarios[AVERAGE].sell_favorable[Wheat] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Wheat] + 150*planning_scenarios[AVERAGE].sell_favorable[Corn] + 0.0*planning_scenarios[AVERAGE].sell_unfavorable[Corn] + 36*planning_scenarios[AVERAGE].sell_favorable[Beats] + 10.0*planning_scenarios[AVERAGE].sell_unfavorable[Beats])
2 Constraint Declarations
crop_min_constraint : Size=3, Index=planning_scenarios[AVERAGE].CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 20.0*planning_scenarios[AVERAGE].plant[Beats] + planning_scenarios[AVERAGE].buy[Beats] - planning_scenarios[AVERAGE].sell_unfavorable[Beats] - planning_scenarios[AVERAGE].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 3.0*planning_scenarios[AVERAGE].plant[Corn] + planning_scenarios[AVERAGE].buy[Corn] - planning_scenarios[AVERAGE].sell_unfavorable[Corn] - planning_scenarios[AVERAGE].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 2.5*planning_scenarios[AVERAGE].plant[Wheat] + planning_scenarios[AVERAGE].buy[Wheat] - planning_scenarios[AVERAGE].sell_unfavorable[Wheat] - planning_scenarios[AVERAGE].sell_favorable[Wheat] : +Inf : True
total_area : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : planning_scenarios[AVERAGE].plant[Wheat] + planning_scenarios[AVERAGE].plant[Corn] + planning_scenarios[AVERAGE].plant[Beats] : 500.0 : True
15 Declarations: CROPS cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible plant buy sell_favorable sell_unfavorable total_area crop_min_constraint net_profit
planning_scenarios[BELOW] : Active=True
1 Set Declarations
CROPS : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'Wheat', 'Corn', 'Beats'}
7 Param Declarations
cost : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 260
Corn : 230
Wheat : 150
crop_yield : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 16.0
Corn : 2.4000000000000004
Wheat : 2.0
max_possible : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 6000.0
Corn : inf
Wheat : inf
min_required : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 240.0
Wheat : 200.0
purchase_price : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 0.0
Corn : 210.0
Wheat : 238.0
sell_price_favorable : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 36
Corn : 150
Wheat : 170
sell_price_unfavorable : Size=3, Index=planning_scenarios[BELOW].CROPS, Domain=NonNegativeReals, Default=None, Mutable=False
Key : Value
Beats : 10.0
Corn : 0.0
Wheat : 0.0
4 Var Declarations
buy : Size=3, Index=planning_scenarios[BELOW].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : 0 : None : True : False : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
plant : Size=3, Index=planning_scenarios[BELOW].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_favorable : Size=3, Index=planning_scenarios[BELOW].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : 6000.0 : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
sell_unfavorable : Size=3, Index=planning_scenarios[BELOW].CROPS
Key : Lower : Value : Upper : Fixed : Stale : Domain
Beats : 0 : None : None : False : True : NonNegativeReals
Corn : 0 : None : None : False : True : NonNegativeReals
Wheat : 0 : None : None : False : True : NonNegativeReals
1 Expression Declarations
net_profit : Size=1, Index=None
Key : Expression
None : - (150*planning_scenarios[BELOW].plant[Wheat] + 230*planning_scenarios[BELOW].plant[Corn] + 260*planning_scenarios[BELOW].plant[Beats]) - (238.0*planning_scenarios[BELOW].buy[Wheat] + 210.0*planning_scenarios[BELOW].buy[Corn] + 0.0*planning_scenarios[BELOW].buy[Beats]) + (170*planning_scenarios[BELOW].sell_favorable[Wheat] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Wheat] + 150*planning_scenarios[BELOW].sell_favorable[Corn] + 0.0*planning_scenarios[BELOW].sell_unfavorable[Corn] + 36*planning_scenarios[BELOW].sell_favorable[Beats] + 10.0*planning_scenarios[BELOW].sell_unfavorable[Beats])
2 Constraint Declarations
crop_min_constraint : Size=3, Index=planning_scenarios[BELOW].CROPS, Active=True
Key : Lower : Body : Upper : Active
Beats : 0.0 : 16.0*planning_scenarios[BELOW].plant[Beats] + planning_scenarios[BELOW].buy[Beats] - planning_scenarios[BELOW].sell_unfavorable[Beats] - planning_scenarios[BELOW].sell_favorable[Beats] : +Inf : True
Corn : 240.0 : 2.4000000000000004*planning_scenarios[BELOW].plant[Corn] + planning_scenarios[BELOW].buy[Corn] - planning_scenarios[BELOW].sell_unfavorable[Corn] - planning_scenarios[BELOW].sell_favorable[Corn] : +Inf : True
Wheat : 200.0 : 2.0*planning_scenarios[BELOW].plant[Wheat] + planning_scenarios[BELOW].buy[Wheat] - planning_scenarios[BELOW].sell_unfavorable[Wheat] - planning_scenarios[BELOW].sell_favorable[Wheat] : +Inf : True
total_area : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : planning_scenarios[BELOW].plant[Wheat] + planning_scenarios[BELOW].plant[Corn] + planning_scenarios[BELOW].plant[Beats] : 500.0 : True
15 Declarations: CROPS cost sell_price_favorable sell_price_unfavorable purchase_price crop_yield min_required max_possible plant buy sell_favorable sell_unfavorable total_area crop_min_constraint net_profit
5 Declarations: planning_scenarios SCENARIOS CROPS enforce_same_planting obj
solver = SolverFactory('ipopt')
solver.solve(m, tee=True)
Ipopt 3.13.2:
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit http://projects.coin-or.org/Ipopt
This version of Ipopt was compiled from source code available at
https://github.com/IDAES/Ipopt as part of the Institute for the Design of
Advanced Energy Systems Process Systems Engineering Framework (IDAES PSE
Framework) Copyright (c) 2018-2019. See https://github.com/IDAES/idaes-pse.
This version of Ipopt was compiled using HSL, a collection of Fortran codes
for large-scale scientific computation. All technical papers, sales and
publicity material resulting from use of the HSL codes within IPOPT must
contain the following acknowledgement:
HSL, a collection of Fortran codes for large-scale scientific
computation. See http://www.hsl.rl.ac.uk.
******************************************************************************
This is Ipopt version 3.13.2, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 12
Number of nonzeros in inequality constraint Jacobian.: 42
Number of nonzeros in Lagrangian Hessian.............: 0
Total number of variables............................: 33
variables with only lower bounds: 30
variables with lower and upper bounds: 3
variables with only upper bounds: 0
Total number of equality constraints.................: 6
Total number of inequality constraints...............: 12
inequality constraints with only lower bounds: 9
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 3
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 7.2199928e+00 2.40e+02 5.11e+01 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 1.8002807e+01 2.40e+02 5.05e+01 -1.0 2.35e+02 - 4.57e-03 1.01e-02h 1
2 4.7679399e+01 2.40e+02 5.09e+01 -1.0 1.58e+02 - 3.92e-04 9.48e-04h 1
3 1.2374072e+02 2.39e+02 5.16e+01 -1.0 1.57e+02 - 6.43e-04 2.12e-03h 1
4 3.9906486e+02 2.37e+02 5.25e+01 -1.0 1.55e+02 - 2.37e-03 7.65e-03h 1
5 5.1604153e+02 2.36e+02 5.23e+01 -1.0 1.63e+02 - 1.61e-03 4.16e-03h 1
6 6.9269623e+02 2.35e+02 5.19e+01 -1.0 1.59e+02 - 1.04e-02 6.17e-03h 1
7 1.1677644e+03 2.30e+02 5.13e+01 -1.0 1.80e+02 - 1.05e-04 2.00e-02h 1
8 1.0569204e+03 2.30e+02 5.13e+01 -1.0 2.74e+04 - 5.09e-05 3.83e-05f 1
9 1.0572289e+03 2.30e+02 5.12e+01 -1.0 2.76e+02 - 6.33e-03 5.82e-05h 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 1.1290620e+03 2.27e+02 5.07e+01 -1.0 3.13e+02 - 6.41e-05 1.46e-02h 1
11 1.1273053e+03 2.27e+02 5.06e+01 -1.0 1.53e+03 - 4.36e-03 2.03e-05f 1
12 1.5097453e+03 2.22e+02 4.99e+01 -1.0 1.90e+02 - 6.58e-03 1.99e-02h 1
13 1.1511177e+03 2.19e+02 4.94e+01 -1.0 5.67e+02 - 6.66e-05 1.46e-02f 1
14 -1.2871873e+05 2.16e+02 4.91e+01 -1.0 1.12e+05 - 6.23e-05 1.27e-02f 1
15 -1.2930322e+05 2.14e+02 4.63e+01 -1.0 8.91e+02 - 1.27e-01 1.06e-02f 1
16 -1.2000705e+05 1.84e+02 4.17e+01 -1.0 4.45e+02 - 3.53e-02 1.38e-01h 1
17 -1.0737267e+05 1.44e+02 3.56e+01 -1.0 3.86e+02 - 3.71e-02 2.18e-01h 1
18 -7.4727388e+04 0.00e+00 7.33e+01 -1.0 1.65e+03 - 2.37e-03 1.00e+00h 1
19 -1.0791818e+05 0.00e+00 6.77e+01 -1.0 2.10e+04 - 7.67e-02 2.05e-01f 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
20 -1.0826496e+05 0.00e+00 6.40e+01 -1.0 3.21e+02 - 5.55e-02 1.34e-01f 1
21 -1.0838277e+05 0.00e+00 6.39e+01 -1.0 1.32e+04 - 8.18e-04 4.58e-03f 1
22 -1.0838396e+05 0.00e+00 2.72e+01 -1.0 5.57e+02 - 5.75e-01 1.09e-03f 1
23 -1.0838720e+05 0.00e+00 1.87e+00 -1.0 6.28e-01 - 9.90e-01 6.64e-01f 1
24 -1.0838710e+05 0.00e+00 1.00e-06 -1.0 4.65e-02 - 1.00e+00 1.00e+00f 1
25 -1.0838992e+05 0.00e+00 2.83e-08 -2.5 1.01e-01 - 1.00e+00 1.00e+00f 1
26 -1.0839000e+05 0.00e+00 1.50e-09 -3.8 2.82e-03 - 1.00e+00 1.00e+00f 1
27 -1.0839000e+05 0.00e+00 1.85e-11 -5.7 1.55e-04 - 1.00e+00 1.00e+00f 1
28 -1.0839000e+05 1.42e-14 4.58e-14 -8.6 1.92e-06 - 1.00e+00 1.00e+00f 1
Number of Iterations....: 28
(scaled) (unscaled)
Objective...............: -1.0839000259619651e+05 -1.0839000259619651e+05
Dual infeasibility......: 4.5776097145183625e-14 4.5776097145183625e-14
Constraint violation....: 1.4210854715202004e-14 1.4210854715202004e-14
Complementarity.........: 2.5067947717597657e-09 2.5067947717597657e-09
Overall NLP error.......: 2.5067947717597657e-09 2.5067947717597657e-09
Number of objective function evaluations = 29
Number of objective gradient evaluations = 29
Number of equality constraint evaluations = 29
Number of inequality constraint evaluations = 29
Number of equality constraint Jacobian evaluations = 29
Number of inequality constraint Jacobian evaluations = 29
Number of Lagrangian Hessian evaluations = 28
Total CPU secs in IPOPT (w/o function evaluations) = 0.003
Total CPU secs in NLP function evaluations = 0.000
EXIT: Optimal Solution Found.
{'Problem': [{'Lower bound': -inf, 'Upper bound': inf, 'Number of objectives': 1, 'Number of constraints': 18, 'Number of variables': 33, 'Sense': 'unknown'}], 'Solver': [{'Status': 'ok', 'Message': 'Ipopt 3.13.2\\x3a Optimal Solution Found', 'Termination condition': 'optimal', 'Id': 0, 'Error rc': 0, 'Time': 0.023143768310546875}], 'Solution': [OrderedDict([('number of solutions', 0), ('number of solutions displayed', 0)])]}
for i, name in enumerate(scenarios.keys()):
print(f"*** Solving scenario {i} = {name} ***")
print_solution(m.planning_scenarios[name])
print("\n")
*** Solving scenario 0 = BELOW ***
Planting decisions:
Wheat: 170.0 acres
Corn: 80.0 acres
Beats: 250.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 48.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 140.0 T
Corn: 0.0 T
Beats: 4000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Solving scenario 1 = AVERAGE ***
Planting decisions:
Wheat: 170.0 acres
Corn: 80.0 acres
Beats: 250.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 225.0 T
Corn: 0.0 T
Beats: 5000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T
*** Solving scenario 2 = ABOVE ***
Planting decisions:
Wheat: 170.0 acres
Corn: 80.0 acres
Beats: 250.0 acres
Purchase decisions:
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0 T
Sales decisions (favorable price):
Wheat: 310.0 T
Corn: 48.0 T
Beats: 6000.0 T
Sales decisions (unfavorable price):
Wheat: 0.0 T
Corn: 0.0 T
Beats: 0.0 T