From b2ed0f67c1e0391657937529d2c92032d33285dc Mon Sep 17 00:00:00 2001 From: Jun He Date: Fri, 31 Mar 2023 15:11:37 -0400 Subject: [PATCH] Added the profiled units --- src/UnitCommitment.jl | 1 + src/instance/migrate.jl | 12 ++ src/instance/read.jl | 185 ++++++++++++++++----------- src/instance/structs.jl | 11 ++ src/model/build.jl | 3 + src/model/formulations/base/punit.jl | 32 +++++ src/solution/solution.jl | 39 ++++-- test/instance/read_test.jl | 22 +++- 8 files changed, 215 insertions(+), 90 deletions(-) create mode 100644 src/model/formulations/base/punit.jl diff --git a/src/UnitCommitment.jl b/src/UnitCommitment.jl index 5d0bf24..3191395 100644 --- a/src/UnitCommitment.jl +++ b/src/UnitCommitment.jl @@ -32,6 +32,7 @@ include("model/formulations/base/psload.jl") include("model/formulations/base/sensitivity.jl") include("model/formulations/base/system.jl") include("model/formulations/base/unit.jl") +include("model/formulations/base/punit.jl") include("model/formulations/CarArr2006/pwlcosts.jl") include("model/formulations/DamKucRajAta2016/ramp.jl") include("model/formulations/Gar1962/pwlcosts.jl") diff --git a/src/instance/migrate.jl b/src/instance/migrate.jl index becb912..317a309 100644 --- a/src/instance/migrate.jl +++ b/src/instance/migrate.jl @@ -17,6 +17,7 @@ function _migrate(json) end version = VersionNumber(version) version >= v"0.3" || _migrate_to_v03(json) + version >= v"0.4" || _migrate_to_v04(json) return end @@ -36,3 +37,14 @@ function _migrate_to_v03(json) end end end + +function _migrate_to_v04(json) + # Migrate thermal units + if json["Generators"] !== nothing + for (gen_name, gen) in json["Generators"] + if gen["Type"] === nothing + gen["Type"] = "Thermal" + end + end + end +end diff --git a/src/instance/read.jl b/src/instance/read.jl index b823dc1..635be8d 100644 --- a/src/instance/read.jl +++ b/src/instance/read.jl @@ -135,6 +135,7 @@ function _from_json(json; repair = true)::UnitCommitmentScenario lines = TransmissionLine[] loads = PriceSensitiveLoad[] reserves = Reserve[] + profiled_units = ProfiledUnit[] function scalar(x; default = nothing) x !== nothing || return default @@ -182,6 +183,7 @@ function _from_json(json; repair = true)::UnitCommitmentScenario timeseries(dict["Load (MW)"]), Unit[], PriceSensitiveLoad[], + ProfiledUnit[], ) name_to_bus[bus_name] = bus push!(buses, bus) @@ -207,90 +209,119 @@ function _from_json(json; repair = true)::UnitCommitmentScenario # Read units for (unit_name, dict) in json["Generators"] + # Read and validate unit type + unit_type = scalar(dict["Type"], default = nothing) + unit_type !== nothing || error("unit $unit_name has no type specified") bus = name_to_bus[dict["Bus"]] - # Read production cost curve - K = length(dict["Production cost curve (MW)"]) - curve_mw = hcat( - [timeseries(dict["Production cost curve (MW)"][k]) for k in 1:K]..., - ) - curve_cost = hcat( - [timeseries(dict["Production cost curve (\$)"][k]) for k in 1:K]..., - ) - min_power = curve_mw[:, 1] - max_power = curve_mw[:, K] - min_power_cost = curve_cost[:, 1] - segments = CostSegment[] - for k in 2:K - amount = curve_mw[:, k] - curve_mw[:, k-1] - cost = (curve_cost[:, k] - curve_cost[:, k-1]) ./ amount - replace!(cost, NaN => 0.0) - push!(segments, CostSegment(amount, cost)) - end - - # Read startup costs - startup_delays = scalar(dict["Startup delays (h)"], default = [1]) - startup_costs = scalar(dict["Startup costs (\$)"], default = [0.0]) - startup_categories = StartupCategory[] - for k in 1:length(startup_delays) - push!( - startup_categories, - StartupCategory( - startup_delays[k] .* time_multiplier, - startup_costs[k], - ), + if lowercase(unit_type) === "thermal" + # Read production cost curve + K = length(dict["Production cost curve (MW)"]) + curve_mw = hcat( + [ + timeseries(dict["Production cost curve (MW)"][k]) for + k in 1:K + ]..., ) - end + curve_cost = hcat( + [ + timeseries(dict["Production cost curve (\$)"][k]) for + k in 1:K + ]..., + ) + min_power = curve_mw[:, 1] + max_power = curve_mw[:, K] + min_power_cost = curve_cost[:, 1] + segments = CostSegment[] + for k in 2:K + amount = curve_mw[:, k] - curve_mw[:, k-1] + cost = (curve_cost[:, k] - curve_cost[:, k-1]) ./ amount + replace!(cost, NaN => 0.0) + push!(segments, CostSegment(amount, cost)) + end - # Read reserve eligibility - unit_reserves = Reserve[] - if "Reserve eligibility" in keys(dict) - unit_reserves = - [name_to_reserve[n] for n in dict["Reserve eligibility"]] - end + # Read startup costs + startup_delays = scalar(dict["Startup delays (h)"], default = [1]) + startup_costs = scalar(dict["Startup costs (\$)"], default = [0.0]) + startup_categories = StartupCategory[] + for k in 1:length(startup_delays) + push!( + startup_categories, + StartupCategory( + startup_delays[k] .* time_multiplier, + startup_costs[k], + ), + ) + end - # Read and validate initial conditions - initial_power = scalar(dict["Initial power (MW)"], default = nothing) - initial_status = scalar(dict["Initial status (h)"], default = nothing) - if initial_power === nothing - initial_status === nothing || - error("unit $unit_name has initial status but no initial power") - else - initial_status !== nothing || - error("unit $unit_name has initial power but no initial status") - initial_status != 0 || - error("unit $unit_name has invalid initial status") - if initial_status < 0 && initial_power > 1e-3 - error("unit $unit_name has invalid initial power") + # Read reserve eligibility + unit_reserves = Reserve[] + if "Reserve eligibility" in keys(dict) + unit_reserves = + [name_to_reserve[n] for n in dict["Reserve eligibility"]] end - initial_status *= time_multiplier - end - unit = Unit( - unit_name, - bus, - max_power, - min_power, - timeseries(dict["Must run?"], default = [false for t in 1:T]), - min_power_cost, - segments, - scalar(dict["Minimum uptime (h)"], default = 1) * time_multiplier, - scalar(dict["Minimum downtime (h)"], default = 1) * time_multiplier, - scalar(dict["Ramp up limit (MW)"], default = 1e6), - scalar(dict["Ramp down limit (MW)"], default = 1e6), - scalar(dict["Startup limit (MW)"], default = 1e6), - scalar(dict["Shutdown limit (MW)"], default = 1e6), - initial_status, - initial_power, - startup_categories, - unit_reserves, - ) - push!(bus.units, unit) - for r in unit_reserves - push!(r.units, unit) + # Read and validate initial conditions + initial_power = + scalar(dict["Initial power (MW)"], default = nothing) + initial_status = + scalar(dict["Initial status (h)"], default = nothing) + if initial_power === nothing + initial_status === nothing || error( + "unit $unit_name has initial status but no initial power", + ) + else + initial_status !== nothing || error( + "unit $unit_name has initial power but no initial status", + ) + initial_status != 0 || + error("unit $unit_name has invalid initial status") + if initial_status < 0 && initial_power > 1e-3 + error("unit $unit_name has invalid initial power") + end + initial_status *= time_multiplier + end + + unit = Unit( + unit_name, + bus, + max_power, + min_power, + timeseries(dict["Must run?"], default = [false for t in 1:T]), + min_power_cost, + segments, + scalar(dict["Minimum uptime (h)"], default = 1) * + time_multiplier, + scalar(dict["Minimum downtime (h)"], default = 1) * + time_multiplier, + scalar(dict["Ramp up limit (MW)"], default = 1e6), + scalar(dict["Ramp down limit (MW)"], default = 1e6), + scalar(dict["Startup limit (MW)"], default = 1e6), + scalar(dict["Shutdown limit (MW)"], default = 1e6), + initial_status, + initial_power, + startup_categories, + unit_reserves, + ) + push!(bus.units, unit) + for r in unit_reserves + push!(r.units, unit) + end + name_to_unit[unit_name] = unit + push!(units, unit) + elseif lowercase(unit_type) === "profiled" + bus = name_to_bus[dict["Bus"]] + pu = ProfiledUnit( + unit_name, + bus, + timeseries(dict["Maximum power (MW)"]), + timeseries(dict["Cost (\$/MW)"]), + ) + push!(bus.profiled_units, pu) + push!(profiled_units, pu) + else + error("unit $unit_name has an invalid type") end - name_to_unit[unit_name] = unit - push!(units, unit) end # Read transmission lines @@ -371,6 +402,8 @@ function _from_json(json; repair = true)::UnitCommitmentScenario time = T, units_by_name = Dict(g.name => g for g in units), units = units, + profiled_units_by_name = Dict(pu.name => pu for pu in profiled_units), + profiled_units = profiled_units, isf = spzeros(Float64, length(lines), length(buses) - 1), lodf = spzeros(Float64, length(lines), length(lines)), ) diff --git a/src/instance/structs.jl b/src/instance/structs.jl index 0c72a23..88a64fd 100644 --- a/src/instance/structs.jl +++ b/src/instance/structs.jl @@ -8,6 +8,7 @@ mutable struct Bus load::Vector{Float64} units::Vector price_sensitive_loads::Vector + profiled_units::Vector end mutable struct CostSegment @@ -73,6 +74,13 @@ mutable struct PriceSensitiveLoad revenue::Vector{Float64} end +mutable struct ProfiledUnit + name::String + bus::Bus + capacity::Vector{Float64} + cost::Vector{Float64} +end + Base.@kwdef mutable struct UnitCommitmentScenario buses_by_name::Dict{AbstractString,Bus} buses::Vector{Bus} @@ -92,6 +100,8 @@ Base.@kwdef mutable struct UnitCommitmentScenario time::Int units_by_name::Dict{AbstractString,Unit} units::Vector{Unit} + profiled_units_by_name::Dict{AbstractString,ProfiledUnit} + profiled_units::Vector{ProfiledUnit} end Base.@kwdef mutable struct UnitCommitmentInstance @@ -108,6 +118,7 @@ function Base.show(io::IO, instance::UnitCommitmentInstance) print(io, "$(length(sc.lines)) lines, ") print(io, "$(length(sc.contingencies)) contingencies, ") print(io, "$(length(sc.price_sensitive_loads)) price sensitive loads, ") + print(io, "$(length(sc.profiled_units)) profiled units, ") print(io, "$(instance.time) time steps") print(io, ")") return diff --git a/src/model/build.jl b/src/model/build.jl index fdfc27e..1bd1987 100644 --- a/src/model/build.jl +++ b/src/model/build.jl @@ -96,6 +96,9 @@ function build_model(; for g in sc.units _add_unit_dispatch!(model, g, formulation, sc) end + for pu in sc.profiled_units + _add_profiled_unit!(model, pu, sc) + end _add_system_wide_eqs!(model, sc) end @objective(model, Min, model[:obj]) diff --git a/src/model/formulations/base/punit.jl b/src/model/formulations/base/punit.jl new file mode 100644 index 0000000..ced3ba9 --- /dev/null +++ b/src/model/formulations/base/punit.jl @@ -0,0 +1,32 @@ +# UnitCommitment.jl: Optimization Package for Security-Constrained Unit Commitment +# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved. +# Released under the modified BSD license. See COPYING.md for more details. + +function _add_profiled_unit!( + model::JuMP.Model, + pu::ProfiledUnit, + sc::UnitCommitmentScenario, +)::Nothing + punits = _init(model, :prod_profiled) + net_injection = _init(model, :expr_net_injection) + for t in 1:model[:instance].time + # Decision variable + punits[sc.name, pu.name, t] = + @variable(model, lower_bound = 0, upper_bound = pu.capacity[t]) + + # Objective function terms + add_to_expression!( + model[:obj], + punits[sc.name, pu.name, t], + pu.cost[t] * sc.probability, + ) + + # Net injection + add_to_expression!( + net_injection[sc.name, pu.bus.name, t], + punits[sc.name, pu.name, t], + 1.0, + ) + end + return +end diff --git a/src/solution/solution.jl b/src/solution/solution.jl index b170f30..25666a2 100644 --- a/src/solution/solution.jl +++ b/src/solution/solution.jl @@ -65,19 +65,22 @@ function solution(model::JuMP.Model)::OrderedDict sol = OrderedDict() for sc in instance.scenarios sol[sc.name] = OrderedDict() - sol[sc.name]["Production (MW)"] = - OrderedDict(g.name => production(g, sc) for g in sc.units) - sol[sc.name]["Production cost (\$)"] = - OrderedDict(g.name => production_cost(g, sc) for g in sc.units) - sol[sc.name]["Startup cost (\$)"] = - OrderedDict(g.name => startup_cost(g, sc) for g in sc.units) - sol[sc.name]["Is on"] = timeseries(model[:is_on], sc.units) - sol[sc.name]["Switch on"] = timeseries(model[:switch_on], sc.units) - sol[sc.name]["Switch off"] = timeseries(model[:switch_off], sc.units) - sol[sc.name]["Net injection (MW)"] = - timeseries(model[:net_injection], sc.buses, sc = sc) - sol[sc.name]["Load curtail (MW)"] = - timeseries(model[:curtail], sc.buses, sc = sc) + if !isempty(sc.units) + sol[sc.name]["Production (MW)"] = + OrderedDict(g.name => production(g, sc) for g in sc.units) + sol[sc.name]["Production cost (\$)"] = + OrderedDict(g.name => production_cost(g, sc) for g in sc.units) + sol[sc.name]["Startup cost (\$)"] = + OrderedDict(g.name => startup_cost(g, sc) for g in sc.units) + sol[sc.name]["Is on"] = timeseries(model[:is_on], sc.units) + sol[sc.name]["Switch on"] = timeseries(model[:switch_on], sc.units) + sol[sc.name]["Switch off"] = + timeseries(model[:switch_off], sc.units) + sol[sc.name]["Net injection (MW)"] = + timeseries(model[:net_injection], sc.buses, sc = sc) + sol[sc.name]["Load curtail (MW)"] = + timeseries(model[:curtail], sc.buses, sc = sc) + end if !isempty(sc.lines) sol[sc.name]["Line overflow (MW)"] = timeseries(model[:overflow], sc.lines, sc = sc) @@ -86,6 +89,16 @@ function solution(model::JuMP.Model)::OrderedDict sol[sc.name]["Price-sensitive loads (MW)"] = timeseries(model[:loads], sc.price_sensitive_loads, sc = sc) end + if !isempty(sc.profiled_units) + sol[sc.name]["Profiled production (MW)"] = + timeseries(model[:prod_profiled], sc.profiled_units, sc = sc) + sol[sc.name]["Profiled production cost (\$)"] = OrderedDict( + pu.name => [ + value(model[:prod_profiled][sc.name, pu.name, t]) * + pu.cost[t] for t in 1:instance.time + ] for pu in sc.profiled_units + ) + end sol[sc.name]["Spinning reserve (MW)"] = OrderedDict( r.name => OrderedDict( g.name => [ diff --git a/test/instance/read_test.jl b/test/instance/read_test.jl index 5fa4a93..81e16d1 100644 --- a/test/instance/read_test.jl +++ b/test/instance/read_test.jl @@ -9,7 +9,7 @@ using UnitCommitment, LinearAlgebra, Cbc, JuMP, JSON, GZip @test repr(instance) == ( "UnitCommitmentInstance(1 scenarios, 6 units, 14 buses, " * - "20 lines, 19 contingencies, 1 price sensitive loads, 4 time steps)" + "20 lines, 19 contingencies, 1 price sensitive loads, 0 profiled units, 4 time steps)" ) @test length(instance.scenarios) == 1 @@ -131,3 +131,23 @@ end @test unit.startup_categories[3].delay == 6 @test unit.initial_status == -200 end + +@testset "read_benchmark profiled-units" begin + instance = UnitCommitment.read("$FIXTURES/case14-profiled.json.gz") + sc = instance.scenarios[1] + @test length(sc.profiled_units) == 2 + + first_pu = sc.profiled_units[1] + @test first_pu.name == "g7" + @test first_pu.bus.name == "b4" + @test first_pu.cost == [100.0 for t in 1:4] + @test first_pu.capacity == [100.0 for t in 1:4] + @test sc.profiled_units_by_name["g7"].name == "g7" + + second_pu = sc.profiled_units[2] + @test second_pu.name == "g8" + @test second_pu.bus.name == "b5" + @test second_pu.cost == [50.0 for t in 1:4] + @test second_pu.capacity == [120.0 for t in 1:4] + @test sc.profiled_units_by_name["g8"].name == "g8" +end