Economic Dispatch

Lecture 15

October 23, 2024

Review and Questions

Last Class

  • Overview of power systems decision problems
  • Many are LPs (at least classically…)
  • Generating capacity expansion

Questions?

Poll Everywhere QR Code

Text: VSRIKRISH to 22333

URL: https://pollev.com/vsrikrish

See Results

Electric Power System Decision Problems

Overview of Electric Power Systems

Power Systems Schematic

Source: Wikipedia

Decisions Problems for Power Systems

Decision Problems for Power Systems by Time Scale

Adapted from Perez-Arriaga, Ignacio J., Hugh Rudnick, and Michel Rivier (2009)

Single-Period Economic Dispatch

Economic Dispatch

Decision Problem: Given a fleet of (online) generators, how do we meet demand at lowest cost?

New Constraints: Power plants are generally subject to engineering constraints that we had previously neglected:

  • Ramping limits
  • Minimum/Maximum power outputs
  • May include network constraints (we will ignore here)

Single-Period Economic Dispatch

What are our variables?

Variable Meaning
\(d\) demand (MWh)
\(y_g\) generation (MWh) by generator \(g \in \mathcal{G}\)
\(VarCost_g\) variable generation cost ($/MWh) for generator \(g\)
\(P^{\text{min/max}}_g\) generation limits (MW) for generator \(g\)

Note on Variable Costs

In practice, variable costs come from:

  • labor costs;
  • equipment upkeep;
  • fuel costs (the big one!).

Fuel is the big variable. The translation of fuel costs to generation costs depends on the efficiency (the heat rate) of the plant.

Note on Variable Costs

These costs are often actually quadratic, not linear, due to efficiency changes.

But we can assume a piecewise linear approximation.

Typical cost curve for a thermal generator

Source: Ross Baldrick, UT Austin

Single-Period Economic Dispatch

Then the economic dispatch problem becomes:

\[ \begin{align} \min_{y_g} & \sum_g VarCost_g \times y_g & \\ \text{subject to:} \quad & \sum_g y_g \geq d & \forall g \in \mathcal{G} \\[0.5em] & y_g \leq P^{\text{max}}_g & \forall g \in \mathcal{G} \\[0.5em] & y_g \geq P^{\text{min}}_g & \forall g \in \mathcal{G} \end{align} \]

Single-Period Example: Data

  • 1 biomass, 50 MW capacity, $5/MWh
  • 1 hydroelectric, 250 MW capacity, $0/MWh
  • 5 natural gas CCGT, 25-220 MW minimum, 50-500 MW capacity $22-36/MWh
  • 7 natural gas CT, 0-73 MW minimum, 48-100 MW capacity, $38-46/MWh

Let’s assume demand is 1600 MWh.

Single-Period Results

Code 
# define sets
# going to ignore renewables for this, so don't look at the last two generators
G = 1:nrow(gens)-2
NSECost = 9000
d = 1600

single_ed = Model(HiGHS.Optimizer)
@variable(single_ed, gens[g, :Pmin] <= y[g in G] <= gens[g, :Pmax])
@objective(single_ed, Min, sum(gens[G, :VarCost] .* y))
@constraint(single_ed, demand, sum(y) >= d)
set_silent(single_ed)
optimize!(single_ed)

leg_gps = repeat(["Dispatched Generation", "Minimum Generation", "Maximum Generation"], inner = nrow(gens)-2)
groupedbar(repeat(gens[G, :Plant], outer=3), Matrix(hcat(Vector(value.(y)), gens[G, :Pmin], gens[G, :Pmax])), group=leg_gps, xrotation=45, bottom_margin=25mm, top_margin=5mm, right_margin=5mm, ylabel="Generation (MWh)", legend=:outerbottom)
plot!(size=(1300, 725))

Cost of Marginal Generation

Code 
gens.cap = gens[:, :Pmax] .- gens[:, :Pmin]
gens_nr = sort!(gens[G, :], [:VarCost, order(:cap, rev=true)])

function plot_supply_curve(supply_curve, demand) 
    rectangle(w, h, x, y) = Shape(x .+ [0,w,w,0], y .+ [0,0,h,h]) 

    p = plot(grid=:false, color_palette=:mk_15, size=(700, 400), titlefontsize=20, top_margin=10mm, right_margin=5mm) 
    plot!(legend=:outerright, legendcolumns=1)
    marg_gen = 0 
    marg_price = 0 
    x = 0 

    plot!(p, rectangle(sum(supply_curve.Pmin), 0.5, x, 0), opacity=.5, label = "minimum", color="black") 
    x = x + sum(supply_curve.Pmin) 
    for i in 1:nrow(supply_curve) 
        if supply_curve[i, :VarCost] == 0 
            plot!(p, rectangle(supply_curve[i, :cap], 0.5,x,0), opacity=.5, label = supply_curve[i, :Plant]) 
        else 
            plot!(p, rectangle(supply_curve[i, :cap], supply_curve[i, :VarCost],x,0), opacity=.5, label = supply_curve[i, :Plant]) 
        end 
        if (x < demand) && (x + supply_curve[i,:cap] > demand) 
            marg_gen = i 
            marg_price = supply_curve[i,:VarCost] 
        end 
        x = x + supply_curve[i,:cap] 
    end 
    vline!([demand],linecolor="black",linewidth=3, linestyle=:dash, label = "demand") 

    title!("Dispatch Stack Supply Curve") 
    xlabel!("Capacity (MW)") 
    ylabel!("Marginal Cost (\$\$/MW)") 

    return (p, marg_price) 
end 
p, marg_price = plot_supply_curve(gens_nr, d)
plot!(p, size=(1300, 500))

Cost of Marginal Generation

Code 
hline!(p, [marg_price],linecolor="blue", linestyle=:dot, linewidth=3, label = "Electricity Price") 
annotate!(p, (400, 20, ("Demand shadow price!", :red, 20)))
quiver!([300], [22], quiver=([0], [15]), color=:red)

Shadow Prices Signs

Note that increasing demand by one MWh decreases cost (minimization problem).

The convention in JuMP is that the shadow price is therefore negative, but this is not universal.

So we can find the price with \[\text{Cost} = |\lambda_\text{demand}|.\]

Dispatch Stack and Merit Order

This supply curve (the dispatch stack) gives the merit order.

Dispatch Stack and Merit Order

What might complicate this simple merit ordering based on variable costs?

Multiple-Period Dispatch

Ramping Constraints

Now, let’s consider multiple time periods.

Not only do we need to meet demand at every time period, but we have additional ramping constraints.

Plants can only increase and decrease their output by so much from time to time, by \(R_g\).

Multi-Period Formulation

\[ \begin{align} \min_{y_{g,t}} & \sum_g VarCost_g \times \sum_t y_{g,t} & \\ \text{subject to:} \quad & \sum_g y_{g,t} \geq d_t & \\[0.3em] & y_{g,t} \leq P^{\text{max}}_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & y_{g,t} \geq P^{\text{min}}_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & \color{red}y_{g,t+1} - y_{g, t} \leq R_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & \color{red}y_{g,t} - y_{g, t+1} \leq R_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \end{align} \]

Multi-Period Generator Data

Ramping constraints can vary strongly by generator type, which, combined with costs, influences whether we view generators as base load or peaking.

  • Nuclear plants generally have a very narrow range in which they can operate;
  • Combustion turbine gas plants can ramp from 0-100% very rapidly.

Multi-Period Generator Data

We’ll make this simple for this problem:

  • Biomass, hydroelectric, CT can ramp from 0-100% each hour.
  • CCGT plants can ramp around 40% of maximum capacity.
  • Renewables: wind (300 MW) and solar (500 MW) can ramp from 0-100% but are limited by their resources each hour, expressed as a capacity factor.

Demand Curve

Code 
NY_demand = DataFrame(CSV.File("data/economic_dispatch/2020_hourly_load_NY.csv"))
rename!(NY_demand, :"Time Stamp" => :Date)
d = NY_demand[:, [:Date, :C]]
rename!(d, :C => :Demand)
n = 289 # pick day
T_period = (n*24+1):((n+1)*24)
d = d[T_period, :]
@df d plot(:Date, :Demand, xlabel="Date", ylabel="Demand (MWh)", label=:false,  xrot=45, bottommargin=24mm)
plot!(size=(1200, 500))
Figure 1: Demand for 2020 in NYISO Zone C

Renewable Variability

Code 
ren_var = DataFrame(CSV.File("data/economic_dispatch/gen_variability.csv"))
p = plot(; xlabel="Hour", ylabel="Capacity Factor")
plot!(p, ren_var[!, :Hour], ren_var[!, :Wind], label="Wind", color=:blue, linewidth=4)
plot!(p, ren_var[!, :Hour], ren_var[!, :Solar], label="Solar", color=:orange, linewidth=4)
ylims!((0, 1))
plot!(size=(1300, 500))
Figure 2: Renewable capacity factors for the given day

Multi-Period Formulation (Renewables)

\[ \begin{align} \min_{y_{g,t}} & \sum_g VarCost_g \times \sum_t y_{g,t} & \\ \text{subject to:} \quad & \sum_g y_{g,t} = d_t & \\[0.3em] & \color{red}y_{g,t} \leq P^{\text{max}}_g \times CF_{g,t} & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & y_{g,t} \geq P^{\text{min}}_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & y_{g,t+1} - y_{g, t} \leq R_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \\[0.3em] & y_{g,t} - y_{g, t+1} \leq R_g & \forall t \in \mathcal{T}, g \in \mathcal{G} \end{align} \]

Multi-Period Results

Code 
G = 1:nrow(gens)
T = 1:nrow(d)

multi_ed = Model(HiGHS.Optimizer)
@variable(multi_ed, gens[g, :Pmin] <= y[g in G, t in T] <= gens[g, :Pmax] * ren_var[t, g+1])
@objective(multi_ed, Min, sum(gens[:, :VarCost] .* [sum(y[g, :]) for g in G]))
@constraint(multi_ed, load[t in T], sum(y[:, t]) == d[t, :Demand])
@constraint(multi_ed, rampup[g in G, t in 1:length(T)-1], y[g, t+1] - y[g, t] <= gens[g, :Ramp])
@constraint(multi_ed, rampdown[g in G, t in 1:length(T)-1], y[g, t] - y[g, t+1] <= gens[g, :Ramp])
set_silent(multi_ed)
optimize!(multi_ed)

multi_ed_noramp = Model(HiGHS.Optimizer)
@variable(multi_ed_noramp, gens[g, :Pmin] <= y2[g in G, t in T] <= gens[g, :Pmax] * ren_var[t, g+1])
@objective(multi_ed_noramp, Min, sum(gens[:, :VarCost] .* [sum(y2[g, :]) for g in G]))
@constraint(multi_ed_noramp, load[t in T], sum(y2[:, t]) == d[t, :Demand])
set_silent(multi_ed_noramp)
optimize!(multi_ed_noramp)


gen = value.(y).data 
p = areaplot(gen', 
    label=permutedims(gens[:, :Plant]), 
    xlabel = "Hour", 
    ylabel ="Generated Electricity (MW)", 
    color_palette=:mk_15,
    grid=:false,
    ylim=(0, 1900),
    top_margin=10mm
)
plot!(legend=:outerright, legendcolumns=1)
plot!(p, size=(1200, 500))

Combining Generator Types

Code 
gen_cbm = DataFrame(gen, :auto)
gen_cbm[!, :Resource] = gens.Resource
gen_cbm_sum = combine(groupby(gen_cbm, [:Resource]), Not([:Resource]) .=> sum)

p = areaplot(Matrix(gen_cbm_sum[:, 2:end])', 
           label=permutedims(gen_cbm_sum[:, :Resource]), 
           xlabel = "Hour", 
           ylabel ="Generated Electricity (MW)", 
           color_palette=:mk_15,
           grid=:false,
           ylim=(0, 1900),
           legend = :outerright
       )
plot!(p, size=(1200, 500))
Figure 3: ED Results grouped by generator type

Impact of Ramping Constraints

Ramping constraints increase the cost by $323.

Code 
p = areaplot(Matrix(gen_cbm_sum[:, 2:end])', 
           label=permutedims(gen_cbm_sum[:, :Resource]), 
           xlabel = "Hour", 
           ylabel ="Generated Electricity (MW)", 
           color_palette=:mk_15,
           grid=:false,
           ylim=(0, 1900),
           legend = :left
       )

gen2 = value.(y2).data
gen2_cbm = DataFrame(gen2, :auto)
gen2_cbm[!, :Resource] = gens.Resource
gen2_cbm_sum = combine(groupby(gen2_cbm, [:Resource]), Not([:Resource]) .=> sum)

p2= groupedbar(Matrix(gen_cbm_sum[:, 2:end] .- gen2_cbm_sum[:, 2:end])', 
    bar_position = :stack,
    label = permutedims(gen_cbm_sum[:, :Resource]),
    xlabel = "Hour", 
    ylabel ="Ramping Constraint\nGeneration Increase (MWh)", 
    color_palette=:mk_15,
    grid=:false,
    legend = :outerright,
    top_margin=10mm,
    right_margin=10mm
)
plot!(p2, layout=(2,1), size=(900, 500))
Figure 4: Impact of ramping constraints

The “Duck Curve”

Ramping and minimum generation play a major role in systems with high levels of renewable penetration.

For example, a prominent feature of grids with large solar generation is the “duck curve” (right).

CAISO Duck Curve

Source: Power Magazine

Key Takeaways

Key Takeaways

  • Economic Dispatch: how to generate power from existing generators at least cost.
  • Engineering constraints (ramping, minimum generation) can weaken simple “merit ordering.”
  • Renewables are limited by resource variability (capacity factor).
  • Ramping constraints play a role in “duck curve” behavior of high-renewable systems.

Upcoming Schedule

Next Classes

Next Week: Air Pollution

Assessments

  • Will start Monday by reviewing Prelim 1.
  • HW4 due 10/31.