3 Choosing an uncertainty budget
Problem 03 | Blending with Uncertain Feed Composition | Level 3: Open-ended challenge
3.1 Problem Statement
A plant blends 100 t/d of a liquid intermediate. Impurity must not exceed 7 wt%. Costs and available flows are deterministic. Impurity has no beneficial role, so only upward deviations matter for the quality upper bound. There are no volume changes, storage, or blending losses.
| Feed | Cost (USD/t) | Available (t/d) | Nominal impurity | Maximum upward deviation |
|---|---|---|---|---|
| A | 100 | 60 | 0.02 | 0.01 |
| B | 80 | 60 | 0.08 | 0.02 |
| C | 55 | 60 | 0.14 | 0.03 |
All fractions are mass fractions. Deviations of 0.01 mean one percentage point, not a 1% relative change.
Your task. Use \(\tilde a_i=a_i+d_i u_i\), \(0\le u_i\le1\), and \(\sum_i u_i\le\Gamma\). Solve the reference case \(\Gamma=1.5\), then evaluate \(\Gamma=0,0.5,1,1.5,2,2.5,3\). Recommend a budget and state what evidence would justify it. There is no unique correct recommendation.
3.1.1 Before you code
Does a budget of 1.5 mean that a constraint has a 50% chance of failure? Why should cost be nondecreasing as the budget grows?
3.1.2 Deliverables
- Define decision variables, units, objective, and constraints. State which data and assumptions are fixed.
- Predict one feature of the solution and explain why it should occur.
- Build and solve a Pyomo model. Check balances and bounds independently.
- Interpret the result and investigate the follow-up question below. For an open-ended challenge, state and defend your chosen policy separately from the reference solution.
3.2 Model Formulation
Let \(x_i\) be feed \(i\) in t/d, with \(0\le x_i\le60\). Define \(c_i\) as cost, \(a_i\) as nominal impurity, and \(d_i\) as maximum upward deviation. Let \(Q=100\) t/d and \(L=0.07\).
The mass balance and cost objective are \[\sum_i x_i=Q,\qquad \min C=\sum_i c_i x_i.\] The impurity constraint compares impurity mass flow with \(LQ\). Multiplying a fixed quality limit by the fixed product rate avoids a ratio.
The worst additional impurity is the support function \[W(x,\Gamma)=\max_u\left\{\sum_i d_i x_i u_i:0\le u_i\le1,\ \sum_i u_i\le\Gamma\right\}.\] Its LP dual introduces \(p\ge0\) and \(q_i\ge0\), both in t/d of impurity: \[p+q_i\ge d_i x_i,\qquad \sum_i a_i x_i+\Gamma p+\sum_i q_i\le LQ.\] Minimizing cost subject to these inequalities is the exact robust counterpart for this uncertainty set. For an independent check, sort \(d_i x_i\) from largest to smallest, add the largest \(\lfloor\Gamma\rfloor\) terms, and add the fractional part times the next term. \(\Gamma=0\) gives nominal protection; \(\Gamma=3\) gives full-box protection. A budget alone does not establish a probability of failure.
3.3 Pyomo Implementation
The code below is the model-building portion of models/blending.py. The level argument selects this problem’s assumptions. Run the complete module from the source bundle to reproduce the solution, including its independent checks:
~/.venvs/optim/bin/python -m models.blending 3"""Blend design: nominal, box-robust, and budget-robust quality."""
import itertools
import pyomo.environ as pyo
from models.common import value, solve
COST = {"A": 100, "B": 80, "C": 55} # USD/t
NOMINAL = {"A": 0.02, "B": 0.08, "C": 0.14} # mass fraction
DELTA = {"A": 0.01, "B": 0.02, "C": 0.03}
CAP = {"A": 60, "B": 60, "C": 60} # t/d
TOTAL = 100 # t/d
LIMIT = 0.07 # maximum impurity mass fraction
def build(level, gamma=1.5):
m = pyo.ConcreteModel(name="Blend quality")
m.I = pyo.Set(initialize=list(COST))
m.x = pyo.Var(m.I, domain=pyo.NonNegativeReals, bounds=lambda m, i: (0, CAP[i]))
m.mass = pyo.Constraint(expr=sum(m.x[i] for i in m.I) == TOTAL)
nominal = sum(NOMINAL[i] * m.x[i] for i in m.I)
if level == 1:
m.quality = pyo.Constraint(expr=nominal <= LIMIT * TOTAL)
elif level == 2:
m.quality = pyo.Constraint(
expr=nominal + sum(DELTA[i] * m.x[i] for i in m.I) <= LIMIT * TOTAL
)
else:
m.p = pyo.Var(domain=pyo.NonNegativeReals)
m.q = pyo.Var(m.I, domain=pyo.NonNegativeReals)
m.support = pyo.Constraint(
m.I, rule=lambda m, i: m.p + m.q[i] >= DELTA[i] * m.x[i]
)
m.quality = pyo.Constraint(
expr=nominal + gamma * m.p + sum(m.q[i] for i in m.I) <= LIMIT * TOTAL
)
m.cost = pyo.Objective(expr=sum(COST[i] * m.x[i] for i in m.I))
m._gamma = gamma
return mThe common solve function calls appsi_highs, checks optimal termination before loading a solution, and checks constraint residuals, bounds, and integer domains. After building the model, use:
from models.common import solve
m = solve(build(3))3.4 Optimal Solution
The reference objective is 8,526.3158 USD/d. This value is optimal for the explicitly stated model and data.
| Feed | Flow (t/d) |
|---|---|
| A | 52.6316 |
| B | 26.3158 |
| C | 21.0526 |
| Quantity | Value |
|---|---|
| Nominal impurity (%) | 6.1053 |
| Full-box impurity (%) | 7.7895 |
| Protected impurity (%) | 7.0000 |
| Gamma | 1.5000 |
3.4.1 Check your solution
Check the sorted-term worst case independently of the dual variables. Confirm the two endpoints match Problems 1 and 2.
The automated residual, bound, and integrality audit passed with a maximum violation of 1.11e-16 in model units. Domain-specific checks passed. Model units differ across equations, so this numerical audit does not replace dimensional analysis.
3.4.2 Uncertainty-budget experiment
| gamma | cost |
|---|---|
| 0.0000 | 8125.0000 |
| 0.5000 | 8333.3333 |
| 1.0000 | 8437.5000 |
| 1.5000 | 8526.3158 |
| 2.0000 | 8608.1081 |
| 2.5000 | 8695.6522 |
| 3.0000 | 8771.4286 |
3.5 Brief Discussion
The protected design costs 8,526.32 USD/d, a premium of 401.32 USD/d (4.94%) over the nominal optimum. Full-box impurity is 7.789%. Protection at a smaller budget must not be described as full-box protection.
3.5.1 Extension to investigate
Propose how historical feed assays could inform the uncertainty set. Discuss correlation and whether the box contains physically implausible combinations.