Decision Models and Linear Programming

---

Lecture 13

October 14, 2026

Overview

• Review of Last Class
• Prescriptive Modeling
• Managing Emissions
• Linear Programs
• The Geometry of an LP
• Key Takeaways and Upcoming Schedule
• References

Review of Last Class

Model Calibration and Validation

Both address why we should trust a model’s results, but they are distinct:

• Calibration: selecting model parameters based on matches to observed data;
• Validation: establishing model credibility by comparing model simulations to independent data (empirical validity) and/or analyzing processes/structure (face validity).

Questions?

Text: VSRIKRISH to 22333

URL: https://pollev.com/vsrikrish

See Results

Prescriptive Modeling

Prescriptive Modeling

If we want to make a decision, we are now in the world of prescriptive modeling.

Recall: Prescriptive modeling is intended to specify an action, policy, or decision.

• Descriptive question: “What happens if I do something (or something happens)?”
• Prescriptive question: “What should I do?”

Decision Models

To make a decision, we need to:

• define decision options (or alternatives);
• provide one or more objective functions to assess performance of different alternatives;
• specify constraints to tell us what alternatives are feasible.

Be Careful Not To Overinterpret Optimization Results

The optimal solution of a model is not an optimal solution of a problem unless the model is a perfect representation of the problem, which it never is.

— Ackoff (1979)

Objectives

Typical objectives can include:

• Minimizing costs (or maximizing profits);
• Minimizing environmental impacts;
• Maximizing some other performance metric.

Objectives and Metrics

Objectives are goals, such as “minimize cost” or “maximize environmental quantity”.

Metrics are functions which measure some relevant quantity, such as a specific cost function.

Many different metrics can be used to specify an objective! The objective function is a combination of an objective and a metric.

Optimal Solutions

Suppose the objective function is to maximize \(Z(\mathbf{x})\).

A solution \(\mathbf{x^*}\) is optimal if \(Z(\mathbf{x^*}) \geq Z(\mathbf{x})\) for all feasible \(\mathbf{x}\).

Note: an optimum is not necessarily unique!

Managing Emissions

Receptors Exposed to Multiple Sources

Concentrations from several stacks add: multiple facilities that individually comply can exceed a standard together and the new combined peak can be hard to predict.

One way to comply is to cut emissions, since concentration is proportional to the emission rate.

Raises the question: with several sources, how much should each one cut?

Stack and Receptor Information

standard = 196                                  # 1-hour SO₂ standard  [µg/m³]
stack_location = [(0, 0), (1.5, 200)]           # (downwind [km], crosswind [m])
emission = [50, 40]                             # uncontrolled emissions  [g/s]
receptor_name = ["A", "B", "C"]
receptor_location = [(2.5, 0), (4.5, 120), (4, 250)]
stack_design = (; H=41, u_wind=2.5)             # last week's stack, Class F

	Downwind	Crosswind	SO\(_2\) emitted
Stack 1	0 km	0 m	50 g/s
Stack 2	1.5 km	200 m	40 g/s
Receptor A	2.5 km	0 m
Receptor B	4.5 km	120 m
Receptor C	4 km	250 m

Problem Parameters and Settings

• Both stacks: effective height 41 m.
• Wind: 2.5 m/s, steady, blowing downwind (left to right on the map).
• Stability: class F, a clear night.
• Limit: the 1-hour SO\(_2\) standard, 75 ppb, about 196 µg/m³ at ground level.

Visualizing the Problem

# Ground-level concentration from both stacks at their current emissions
function combined_ugm3(x, y)
    total = 0.0
    for i in 1:length(emission)
        sx, sy = stack_location[i]
        total += plume_ugm3(x - sx, y - sy, 0.0, (; stack_design..., Q=emission[i]))
    end
    return total
end

map_x = -0.3:0.02:6
map_y = -400:5:700
map_C = [combined_ugm3(x, y) for y in map_y, x in map_x]

p_map = contourf(map_x, map_y, map_C, color=:viridis, linewidth=0, clims=(0, 1200),
    colorbar_title="\nC  [µg/m³]", colorbar_titlefontsize=18,
    xlabel="Distance downwind  [km]", ylabel="Crosswind  [m]")
contour!(p_map, map_x, map_y, map_C, levels=[standard], color=:black, linewidth=3,
    colorbar_entry=false)
for i in 1:length(emission)
    sx, sy = stack_location[i]
    scatter!(p_map, [sx], [sy], markershape=:utriangle, markersize=12, color=:white,
        markerstrokecolor=:black)
    annotate!(p_map, sx + 0.35, sy + 110, text("Stack $(i)", 16, :white, plot_font))
end
for j in 1:length(receptor_name)
    rx, ry = receptor_location[j]
    scatter!(p_map, [rx], [ry], markershape=:circle, markersize=11, color=cb_orange,
        markerstrokecolor=:black)
    annotate!(p_map, rx + 0.25, ry, text(receptor_name[j], 20, :white, plot_font))
end
plot!(p_map, size=(1050, 330), left_margin=12mm, right_margin=22mm,
    bottom_margin=12mm, top_margin=5mm)

Figure 1: Ground-level SO₂ from the two stacks with no removal, on a class F night with a 2.5 m/s wind blowing left to right. Black contour: the 1-hour standard.

Goal: meet the 196 µg/m³ standard at all three receptors, at minimum cost.

Decision Variables

What are the decision variables?

\(R_1\), \(R_2\): fraction of SO\(_2\) removed at each stack.

Are there limits to what values they can take?

Each is bounded below by 0 and above by 1:

\[0 \leq R_1 \leq 1, 0 \leq R_2 \leq 1.\]

Objective Function

Removal costs $20,000 per year for each g/s removed at stack 1, and $30,000 at stack 2 (an older plant).

\[\min_{R_1, R_2} \quad 1000 R_1 + 1200 R_2 \qquad [\$1000/\text{yr}]\]

Developing Constraints

We can write the concentration at receptor \(j\) is

\[C_j = t_{1j}\,Q_1 (1 - R_1) + t_{2j}\,Q_2 (1 - R_2),\]

where \(t_{ij}\) is the transmission factor: the concentration at receptor \(j\) per g/s emitted from stack \(i\).

What Is a Transmission Factor?

At ground level (\(z = 0\)), the plume equation is the emission rate \(Q\) times a factor that does not depend on \(Q\):

\[C(x, y) = Q \times \underbrace{\frac{1}{\pi u \sigma_y \sigma_z} \exp\left(-\frac{y^2}{2\sigma_y^2}\right) \exp\left(-\frac{H^2}{2\sigma_z^2}\right)}_{\text{concentration per g/s emitted}}\]

\(t_{ij}\) is that factor, with \(x\) and \(y\) the downwind and crosswind distances from stack \(i\) to receptor \(j\), in µg/m³ per g/s.

Setting z = 0 makes the two exponentials in the plume equation equal (the plume and its reflection off the ground), and their sum, 2 exp(-H^2 / 2 sigma_z^2), cancels the 2 in 2 pi. The sigmas are evaluated at the downwind distance x from that stack. Because the factor does not involve Q, each receptor’s concentration is a weighted sum of the stacks’ emissions, which is what makes this problem linear. A different night (another stability class or wind speed) gives different factors, and so a different problem.

Example: Stack 1 to Receptor A

# Transmission factor from stack 1 to receptor A, one piece at a time
x_1A = receptor_location[1][1] - stack_location[1][1]     # [km]
sy_1A, sz_1A = sig_y(x_1A), sig_z(x_1A)                    # [m]
prefactor_1A = 1e6 / (π * stack_design.u_wind * sy_1A * sz_1A)
height_term_1A = exp(-stack_design.H^2 / (2sz_1A^2))
t_1A = round(prefactor_1A * height_term_1A; digits=1)       # [µg/m³ per g/s], to the 0.1 shown

Receptor A is 2.5 km directly downwind of stack 1 (\(y = 0\)), where \(\sigma_y\) = 77.1 m and \(\sigma_z\) = 25.2 m on a class F night:

\[\begin{aligned} t_{1A} &= \frac{10^6}{\pi \times 2.5 \times 77.1 \times 25.2} \exp\left(-\frac{41^2}{2 \times 25.2^2}\right) \\ &= 65.4 \times 0.267 = 17.5\ \text{µg/m}^3 \text{ per g/s} \end{aligned}\]

So stack 1’s 50 g/s put 50 × 17.5 = 875 µg/m³ at A.

Board item if there is time: work t_1A from the formula two slides back, then show this slide as the check. The height term, about 0.27, is how much lower the ground-level concentration is than if the same plume had been released at ground level: at 2.5 km on a stable night much of it is still aloft.

Transmission Factors

# Concentration at a receptor per g/s emitted from a stack  [µg/m³ per g/s]
function transmission(stack, receptor)
    sx, sy = stack_location[stack]
    rx, ry = receptor_location[receptor]
    return plume_ugm3(rx - sx, ry - sy, 0.0, (; stack_design..., Q=1))
end

n_stacks, n_receptors = length(emission), length(receptor_name)
# Rounded to the 0.1 shown in the table, so that every number that follows
# can be checked by hand from the table
t = [round(transmission(i, j); digits=1) for j in 1:n_receptors, i in 1:n_stacks]   # rows: receptors
uncontrolled = t * emission                                       # [µg/m³]

println("| Receptor | \$t_{1j}\$ (µg/m³ per g/s) | \$t_{2j}\$ (µg/m³ per g/s) | Uncontrolled \$C_j\$ (µg/m³) |")
println("|:--|--:|--:|--:|")
for j in 1:n_receptors
    println("| $(receptor_name[j]) | $(round(t[j, 1]; digits=1)) | $(round(t[j, 2]; digits=1)) | $(Int(round(uncontrolled[j]))) |")
end

Receptor	\(t_{1j}\) (µg/m³ per g/s)	\(t_{2j}\) (µg/m³ per g/s)	Uncontrolled \(C_j\) (µg/m³)
A	17.5	0.0	875
B	9.0	11.5	910
C	1.5	14.2	643

All three receptors are above the standard before any removal.

Each Receptor’s SO\(_2\), by Stack

# Each stack's contribution at each receptor  [µg/m³]
from_stack = zeros(n_receptors, n_stacks)
for j in 1:n_receptors
    for i in 1:n_stacks
        from_stack[j, i] = t[j, i] * emission[i]
    end
end

# Stacked bars: each receptor's total, with stack 1's part drawn over its base
p_sources = bar(1:n_receptors, uncontrolled, color=cb_orange, bar_width=0.55,
    label="from stack 2")
bar!(p_sources, 1:n_receptors, from_stack[:, 1], color=cb_blue, bar_width=0.55,
    label="from stack 1")
hline!(p_sources, [standard], color=:black, linestyle=:dash, linewidth=2.5,
    label="standard")
for j in 1:n_receptors
    annotate!(p_sources, j, uncontrolled[j] + 60,
        text(string(Int(round(uncontrolled[j]))), 16, :black, plot_font))
end
plot!(p_sources, xticks=(1:n_receptors, "Receptor " .* receptor_name),
    ylabel="C  [µg/m³]", ylims=(0, 1050), legend=:outerright,
    size=(1000, 300), left_margin=12mm, bottom_margin=12mm)

Figure 2: SO₂ at each receptor with no removal, split by the stack it comes from. Dashed: the 1-hour standard, which each receptor must be brought down to.

A gets all of its SO\(_2\) from stack 1, so its constraint involves \(R_1\) alone. B gets 49% from stack 1, and C gets 88% from stack 2.

Problem Map (Revisited)

Figure 3: Map of the emissions problem.

The Constraint at Receptor B

Write the constraint that keeps receptor B at or below 196 µg/m³, in terms of \(R_1\) and \(R_2\).

\[450 R_1 + 460 R_2 \geq 714\]

Board item. Start from C_B = t_1B Q_1 (1 - R_1) + t_2B Q_2 (1 - R_2) <= 196. Expand: t_1B Q_1 + t_2B Q_2 - t_1B Q_1 R_1 - t_2B Q_2 R_2 <= 196. The first two terms are the uncontrolled concentration at B. Move the removals to the right and the standard to the left: t_1B Q_1 R_1 + t_2B Q_2 R_2 >= uncontrolled_B - 196. The trap is the direction of the inequality: concentration has an upper limit, so removal has a lower limit.

The Full Problem

\[\begin{alignedat}{3} & \min_{R_1, R_2} &\quad 1000 R_1 + 1200 R_2 & \\ & \text{subject to:} & 875 R_1 &\geq 679 & \quad \text{(C1)} \\ & & 450 R_1 + 460 R_2 &\geq 714 & \quad \text{(C2)} \\ & & 75 R_1 + 568 R_2 &\geq 447 & \quad \text{(C3)} \\ & & 0 \leq R_1 \leq 1, \quad 0 \leq R_2 &\leq 1 & \end{alignedat}\]

C1, C2 and C3 are the limits at receptors A, B and C.

Linear Programs

What Makes a Linear Program?

A linear program (LP), or a linear optimization model, has the following characteristics:

• Linearity: The objective function and constraints are all linear;
• Divisibility: The decision variables are continuous (they can be fractional levels);
• Certainty: The problem is deterministic.

Ask whether the emissions problem qualifies. It does: concentrations are linear in the removals because plumes add, removals are continuous, and the plume model is deterministic. Each of the three properties is also an assumption we could relax later: integer decisions (MIP), and uncertainty (stochastic programming).

Standard Form

Every LP can be written the same way:

\[\max_{\mathbf{x}} \quad \mathbf{c}^\top \mathbf{x} \qquad \text{subject to} \quad \mathbf{A}\mathbf{x} \leq \mathbf{b}, \quad \mathbf{x} \geq \mathbf{0}\]

• A minimization becomes maximizing the negative: \(\min f = -\max(-f)\).
• A \(\geq\) constraint becomes \(\leq\) by multiplying both sides by \(-1\).
• Upper bounds become ordinary constraints.

Solvers and the theory behind them, including Monday’s simplex walk, work from this form.

What Are \(\mathbf{c}\), \(\mathbf{A}\), and \(\mathbf{b}\)?

\[\max_{\mathbf{x}} \quad \mathbf{c}^\top \mathbf{x} \qquad \text{subject to} \quad \mathbf{A}\mathbf{x} \leq \mathbf{b}, \quad \mathbf{x} \geq \mathbf{0}\]

• \(\mathbf{x}\): the decision variables, stacked in a vector.
• \(\mathbf{c}\): the objective coefficients, one for each variable.
• \(\mathbf{A}\): the constraint coefficients, with a row for each constraint and a column for each variable.
• \(\mathbf{b}\): the right-hand sides, one constant for each constraint.

The compact form is only useful once each symbol has a meaning, and anyone working from the slides alone has no board to fall back on. Count before converting: n decision variables give c n entries and A n columns, and each constraint gives one row of A and one entry of b. Every term with a decision variable goes on the left of its constraint, so b holds constants only. Non-negativity is stated separately as x >= 0 and is not a row of A; upper bounds are rows.

Our Problem in Standard Form

Convert the emissions problem to standard form.

\[\begin{alignedat}{3} & \max_{R_1, R_2} &\quad -1000 R_1 - 1200 R_2 & \\ & \text{subject to:} & -875 R_1 &\leq -679 & \quad \text{(C1)} \\ & & -450 R_1 - 460 R_2 &\leq -714 & \quad \text{(C2)} \\ & & -75 R_1 - 568 R_2 &\leq -447 & \quad \text{(C3)} \\ & & R_1 \leq 1, \quad R_2 &\leq 1 & \\ & & R_1, R_2 &\geq 0 & \end{alignedat}\]

Board item. Three moves: negate the objective, multiply each receptor constraint by minus one, and keep the bounds as constraints. Nothing about the solution changes; the optimal removals are the same and the optimal value just changes sign. This is the form Monday’s simplex walk starts from.

Our \(\mathbf{c}\), \(\mathbf{A}\), and \(\mathbf{b}\)

With \(\mathbf{x} = (R_1, R_2)^\top\), each row of \(\mathbf{A}\) and \(\mathbf{b}\) is one constraint:

\[\mathbf{c} = \begin{pmatrix} -1000 \\ -1200 \end{pmatrix} \;\; \mathbf{A} = \begin{pmatrix} -875 & 0 \\ -450 & -460 \\ -75 & -568 \\ 1 & 0 \\ 0 & 1 \end{pmatrix} \;\; \mathbf{b} = \begin{pmatrix} -679 \\ -714 \\ -447 \\ 1 \\ 1 \end{pmatrix} \quad \begin{matrix} \text{(C1)} \\ \text{(C2)} \\ \text{(C3)} \\ R_1 \leq 1 \\ R_2 \leq 1 \end{matrix}\]

\(R_1, R_2 \geq 0\) is already \(\mathbf{x} \geq \mathbf{0}\), so it adds no rows.

Have students read c, A and b off the previous slide before showing this; it can be done at the board instead. The rows follow the order of the constraints on the previous slide. The 0 in the first row is R_2’s coefficient in C1, the constraint at receptor A: a variable that does not appear in a constraint still needs an entry. The two upper bounds were written on one line but are two constraints, so two rows. Row order does not matter as long as A and b use the same one. Homework 6 asks for c, A and b explicitly for its own LP.

The Geometry of an LP

Where Can An Optimum Be?

# A made-up LP: four constraint lines, each x2 = slope * x1 + intercept
sketch_slope = [4.5, -1, -1.5, 0.5]
sketch_intercept = [0, 16, 12, 0]

# Where constraint lines i and k cross
function sketch_crossing(i, k)
    x1 = (sketch_intercept[k] - sketch_intercept[i]) / (sketch_slope[i] - sketch_slope[k])
    return (x1, sketch_slope[i] * x1 + sketch_intercept[i])
end

# The feasible region's corners, going around it
sketch_corners = [sketch_crossing(1, 3), sketch_crossing(1, 2), sketch_crossing(2, 4),
    sketch_crossing(4, 3)]

p_sketch = plot(Shape(first.(sketch_corners), last.(sketch_corners)), color=cb_blue,
    alpha=0.25, linewidth=0, xlims=(0, 14), ylims=(0, 16), ticks=false,
    xlabel=L"x_1", ylabel=L"x_2")
sketch_x1 = 0:0.5:14
for i in 1:length(sketch_slope)
    plot!(p_sketch, sketch_x1, sketch_slope[i] .* sketch_x1 .+ sketch_intercept[i],
        color=:gray30, linewidth=2.5)
end
plot!(p_sketch, size=(580, 480), left_margin=5mm, bottom_margin=5mm)

Figure 4: A sketch of a linear program with two decision variables: each line is a constraint, and the shaded region satisfies all of them.

• Can the optimum be in the interior?
• What about along an edge, but not at a corner?

Board item: redraw the sketch and work on it there, while students annotate their copy. The figure is left bare on purpose: no objective, no level lines, no labels. Pick an objective and draw two or three of its level lines; they are parallel because the objective is linear. From any interior point you can move across the level lines to a better point, so the optimum is on the boundary. Along an edge the objective changes at a constant rate, so one end is better, unless the edge is parallel to the level lines, in which case the whole edge ties and its corners are optimal too. Either way there is an optimum at a corner, which is what the simplex method uses on Monday. The next slides do the same for the emissions problem.

The Feasible Region

# Is a removal plan within the bounds and below the standard everywhere?
is_feasible(r) = all(from_stack * r .>= required .- 1e-6) &&
    all(-1e-9 .<= r .<= 1 + 1e-9)

# Every boundary as (a1, a2, rhs) for a1 R1 + a2 R2 = rhs:
# the receptor constraints C1 to C3, then the four bounds
boundaries = [[(from_stack[j, 1], from_stack[j, 2], required[j]) for j in 1:n_receptors];
              [(1, 0, 0), (1, 0, 1), (0, 1, 0), (0, 1, 1)]]

function crossing(l1, l2)
    determinant = l1[1] * l2[2] - l1[2] * l2[1]
    if abs(determinant) < 1e-9
        return [NaN, NaN]                       # parallel boundaries never meet
    end
    return [(l1[3] * l2[2] - l1[2] * l2[3]) / determinant,
            (l1[1] * l2[3] - l1[3] * l2[1]) / determinant]
end

# Corners: crossings of two boundaries that satisfy every constraint
n_boundaries = length(boundaries)
crossings = [crossing(boundaries[i], boundaries[k]) for i in 1:n_boundaries for k in (i+1):n_boundaries]
corners = unique(p -> round.(p; digits=6), filter(p -> !any(isnan, p) && is_feasible(p), crossings))
centre = sum(corners) / length(corners)
corners = sort(corners, by=p -> atan(p[2] - centre[2], p[1] - centre[1]))

line_style = [:solid, :dash, :dashdot]
function draw_region(r1_max=1.04)
    p = plot(Shape(first.(corners), last.(corners)), color=cb_blue, alpha=0.25,
        linewidth=0, xlims=(0.5, r1_max), ylims=(0.375, 1.05), xlabel=L"R_1",
        ylabel=L"R_2")
    for j in 1:n_receptors
        a1, a2, rhs = boundaries[j]
        # Label each line at a point where it is inside the axes
        if abs(a2) < 1e-3
            vline!(p, [rhs / a1], color=:gray30, linestyle=line_style[j], linewidth=2.5)
            label_at = (rhs / a1 - 0.028, 0.425)
        else
            r1 = 0.5:0.01:r1_max
            plot!(p, r1, (rhs .- a1 .* r1) ./ a2, color=:gray30, linestyle=line_style[j],
                linewidth=2.5)
            r1_at_09 = (rhs - a2 * 0.9) / a1
            label_at = 0.54 <= r1_at_09 <= 0.96 ? (r1_at_09 + 0.048, 0.915) : (0.54, (rhs - a1 * 0.54) / a2 + 0.065)
        end
        annotate!(p, label_at..., text("C$(j)", 20, :gray20, plot_font))
    end
    plot!(p, [0, 1, 1], [1, 1, 0], color=:black, linewidth=2)
    return p
end

p_region = draw_region()
plot!(p_region, size=(900, 380), left_margin=12mm, bottom_margin=12mm)

Figure 5: Removals that meet the standard at all three receptors (shaded).

Each constraint cuts away plans that break one receptor’s limit.

Where Is the Cheapest Plan?

cost = [20, 30] .* emission                      # [$1000/yr per unit fraction removed]
corner_cost = [sum(cost .* c) for c in corners]
optimal_cost, best = findmin(corner_cost)
optimum = corners[best]

# Widen the plot a little so each level line can be labelled at its right end
p_level = draw_region(1.16)
r1 = 0.5:0.01:1.16
levels = [optimal_cost, optimal_cost + 200, optimal_cost + 400]
for k in 1:length(levels)
    plot!(p_level, r1, (levels[k] .- cost[1] .* r1) ./ cost[2], color=cb_vermillion,
        linestyle=:dot, linewidth=2)
    label_y = (levels[k] - cost[1] * 1.104) / cost[2] + 0.04
    annotate!(p_level, 1.104, label_y,
        text("\$$(round(levels[k] / 1000; digits=1))M/yr", 14, :black, plot_font))
end
scatter!(p_level, [optimum[1]], [optimum[2]], markershape=:star5, markersize=16,
    color=cb_vermillion, markerstrokecolor=:black)
annotate!(p_level, optimum[1] - 0.05, optimum[2] - 0.1125,
    text("cheapest plan", 18, :black, plot_font))
plot!(p_level, ylims=(0.375, 1.075), size=(900, 320), left_margin=12mm, bottom_margin=12mm)

Figure 6: Cost level lines (dotted), labelled with their annual cost, slide toward lower cost until they last touch the region, at a corner.

Cheapest, at $1.7M/yr: \(R_1\) = 0.90 (45.0 of stack 1’s 50 g/s removed) and \(R_2\) = 0.67 (26.8 of stack 2’s 40 g/s).

Cost is constant along each dotted line, and the lines are parallel because the cost is linear. Lower cost is toward the origin. Slide the line toward lower cost until it is about to leave the shaded region: the last point it touches is the optimum, and it is a corner. Which corner wins depends on the cost ratio: with these costs, this corner is optimal as long as stack 2’s removals cost between about 1.3 and 9 times stack 1’s. Below that the corner where C1 and C2 meet wins instead. That sensitivity is Monday’s topic.

Reading the Solution

Receptor	Concentration at the cheapest plan
A	84 µg/m³
B	196 µg/m³
C	196 µg/m³

C2 and C3 are binding constraints: receptors B and C sit exactly at the standard. C1 is a slack constraint: the optimal solution is not constrained by it.

When the Cost Isn’t Linear

# Rising marginal cost: each extra g/s removed costs more than the last
curved_cost(r1, r2) = cost[1] * r1^2 + cost[2] * r2^2

# Search a fine grid of plans for the cheapest one that meets the standard
grid_r1 = 0:0.001:1
grid_r2 = 0:0.00125:1
grid_cost = [is_feasible([r1, r2]) ? curved_cost(r1, r2) : Inf for r2 in grid_r2, r1 in grid_r1]
curved_best = argmin(grid_cost)
curved_optimum = [grid_r1[curved_best[2]], grid_r2[curved_best[1]]]

p_curved = draw_region()
contour!(p_curved, 0.5:0.005:1.04, 0.375:0.00625:1.05, (r1, r2) -> curved_cost(r1, r2),
    levels=curved_cost(curved_optimum...) .+ [0, 150, 300], color=cb_vermillion,
    linestyle=:dot, linewidth=2, colorbar_entry=false)
scatter!(p_curved, [curved_optimum[1]], [curved_optimum[2]], markershape=:star5,
    markersize=16, color=cb_vermillion, markerstrokecolor=:black)
annotate!(p_curved, curved_optimum[1] + 0.06, curved_optimum[2] - 0.125,
    text("on an edge", 18, :black, plot_font))
plot!(p_curved, size=(900, 320), left_margin=12mm, bottom_margin=12mm)

Figure 7: With rising marginal costs, cost contours curve and the cheapest plan sits on an edge.

With cost \(1000 R_1^2 + 1200 R_2^2\) ($1000/yr), each g/s costs more than the last, and the optimum leaves the corner.

The curved cost equals the linear one at full removal (1000 x 1 = 1000 x 1^2), so the two are comparable; only the shape differs. Its contours bend, and one of them touches the region partway along C2’s edge rather than at a corner.

Key Takeaways and Upcoming Schedule

Key Takeaways

• A decision model: decision variables (units, bounds), an objective, constraints.
• Constraints often come from fate and transport models or economic/budget constraints.
• An LP will always have an optimum at a corner.
• Every LP can be written in standard form.
• The optimal plan is only as good as the model behind it.

Next Classes

Monday: How to solve an LP and shadow prices.

Wednesday: Quiz 3, then the project proposal peer review.

Assessments

Homework 5: Assigned today, due Thursday, October 22.

Project Proposal: Due Friday, October 16.

Mini-Project 1: Due Friday, October 23.

References

References

Ackoff, R. L. (1979). The Future of Operational Research is Past. J. Oper. Res. Soc., 30, 93–104. https://doi.org/10.1057/jors.1979.22