Added the profiled units

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")

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

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
+        if lowercase(unit_type) === "thermal"
-        K = length(dict["Production cost curve (MW)"])
+            # Read production cost curve
-        curve_mw = hcat(
+            K = length(dict["Production cost curve (MW)"])
-            [timeseries(dict["Production cost curve (MW)"][k]) for k in 1:K]...,
+            curve_mw = hcat(
-        )
+                [
-        curve_cost = hcat(
+                    timeseries(dict["Production cost curve (MW)"][k]) for
-            [timeseries(dict["Production cost curve (\$)"][k]) for k in 1:K]...,
+                    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],
-                ),
             )
-        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
+            # Read startup costs
-        unit_reserves = Reserve[]
+            startup_delays = scalar(dict["Startup delays (h)"], default = [1])
-        if "Reserve eligibility" in keys(dict)
+            startup_costs = scalar(dict["Startup costs (\$)"], default = [0.0])
-            unit_reserves =
+            startup_categories = StartupCategory[]
-                [name_to_reserve[n] for n in dict["Reserve eligibility"]]
+            for k in 1:length(startup_delays)
-        end
+                push!(
+                    startup_categories,
+                    StartupCategory(
+                        startup_delays[k] .* time_multiplier,
+                        startup_costs[k],
+                    ),
+                )
+            end
 
-        # Read and validate initial conditions
+            # Read reserve eligibility
-        initial_power = scalar(dict["Initial power (MW)"], default = nothing)
+            unit_reserves = Reserve[]
-        initial_status = scalar(dict["Initial status (h)"], default = nothing)
+            if "Reserve eligibility" in keys(dict)
-        if initial_power === nothing
+                unit_reserves =
-            initial_status === nothing ||
+                    [name_to_reserve[n] for n in dict["Reserve eligibility"]]
-                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(
+            # Read and validate initial conditions
-            unit_name,
+            initial_power =
-            bus,
+                scalar(dict["Initial power (MW)"], default = nothing)
-            max_power,
+            initial_status =
-            min_power,
+                scalar(dict["Initial status (h)"], default = nothing)
-            timeseries(dict["Must run?"], default = [false for t in 1:T]),
+            if initial_power === nothing
-            min_power_cost,
+                initial_status === nothing || error(
-            segments,
+                    "unit $unit_name has initial status but no initial power",
-            scalar(dict["Minimum uptime (h)"], default = 1) * time_multiplier,
+                )
-            scalar(dict["Minimum downtime (h)"], default = 1) * time_multiplier,
+            else
-            scalar(dict["Ramp up limit (MW)"], default = 1e6),
+                initial_status !== nothing || error(
-            scalar(dict["Ramp down limit (MW)"], default = 1e6),
+                    "unit $unit_name has initial power but no initial status",
-            scalar(dict["Startup limit (MW)"], default = 1e6),
+                )
-            scalar(dict["Shutdown limit (MW)"], default = 1e6),
+                initial_status != 0 ||
-            initial_status,
+                    error("unit $unit_name has invalid initial status")
-            initial_power,
+                if initial_status < 0 && initial_power > 1e-3
-            startup_categories,
+                    error("unit $unit_name has invalid initial power")
-            unit_reserves,
+                end
-        )
+                initial_status *= time_multiplier
-        push!(bus.units, unit)
+            end
-        for r in unit_reserves
+
-            push!(r.units, unit)
+            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)),
     )

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

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])

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

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)"] =
+        if !isempty(sc.units)
-            OrderedDict(g.name => production(g, sc) for g in sc.units)
+            sol[sc.name]["Production (MW)"] =
-        sol[sc.name]["Production cost (\$)"] =
+                OrderedDict(g.name => production(g, sc) for g in sc.units)
-            OrderedDict(g.name => production_cost(g, sc) for g in sc.units)
+            sol[sc.name]["Production cost (\$)"] =
-        sol[sc.name]["Startup cost (\$)"] =
+                OrderedDict(g.name => production_cost(g, sc) for g in sc.units)
-            OrderedDict(g.name => startup_cost(g, sc) for g in sc.units)
+            sol[sc.name]["Startup cost (\$)"] =
-        sol[sc.name]["Is on"] = timeseries(model[:is_on], sc.units)
+                OrderedDict(g.name => startup_cost(g, sc) for g in sc.units)
-        sol[sc.name]["Switch on"] = timeseries(model[:switch_on], sc.units)
+            sol[sc.name]["Is on"] = timeseries(model[:is_on], sc.units)
-        sol[sc.name]["Switch off"] = timeseries(model[:switch_off], sc.units)
+            sol[sc.name]["Switch on"] = timeseries(model[:switch_on], sc.units)
-        sol[sc.name]["Net injection (MW)"] =
+            sol[sc.name]["Switch off"] =
-            timeseries(model[:net_injection], sc.buses, sc = sc)
+                timeseries(model[:switch_off], sc.units)
-        sol[sc.name]["Load curtail (MW)"] =
+            sol[sc.name]["Net injection (MW)"] =
-            timeseries(model[:curtail], sc.buses, sc = sc)
+                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 => [

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