From 7c907a6eb54a3ba4c276a35c0f9f4f406d2b5cfd Mon Sep 17 00:00:00 2001 From: Alinson S Xavier Date: Thu, 29 Jul 2021 11:14:12 -0500 Subject: [PATCH] Implement randomization method from XavQiuAhm2021 --- Project.toml | 3 +- src/UnitCommitment.jl | 2 +- src/transform/randomize.jl | 53 ----- src/transform/randomize/XavQiuAhm2021.jl | 209 ++++++++++++++++++ test/runtests.jl | 4 +- .../transform/randomize/XavQiuAhm2021_test.jl | 63 ++++++ test/transform/randomize_test.jl | 43 ---- 7 files changed, 278 insertions(+), 99 deletions(-) delete mode 100644 src/transform/randomize.jl create mode 100644 src/transform/randomize/XavQiuAhm2021.jl create mode 100644 test/transform/randomize/XavQiuAhm2021_test.jl delete mode 100644 test/transform/randomize_test.jl diff --git a/Project.toml b/Project.toml index e846a41..c7b62f5 100644 --- a/Project.toml +++ b/Project.toml @@ -15,6 +15,7 @@ Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" PackageCompiler = "9b87118b-4619-50d2-8e1e-99f35a4d4d9d" Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7" +Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [compat] @@ -30,8 +31,8 @@ julia = "1" [extras] Cbc = "9961bab8-2fa3-5c5a-9d89-47fab24efd76" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" Gurobi = "2e9cd046-0924-5485-92f1-d5272153d98b" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] test = ["Cbc", "Test", "Gurobi"] diff --git a/src/UnitCommitment.jl b/src/UnitCommitment.jl index 5b727a6..6eb5638 100644 --- a/src/UnitCommitment.jl +++ b/src/UnitCommitment.jl @@ -48,7 +48,7 @@ include("solution/warmstart.jl") include("solution/write.jl") include("transform/initcond.jl") include("transform/slice.jl") -include("transform/randomize.jl") +include("transform/randomize/XavQiuAhm2021.jl") include("utils/log.jl") include("validation/repair.jl") include("validation/validate.jl") diff --git a/src/transform/randomize.jl b/src/transform/randomize.jl deleted file mode 100644 index f8d5bbe..0000000 --- a/src/transform/randomize.jl +++ /dev/null @@ -1,53 +0,0 @@ -# UnitCommitment.jl: Optimization Package for Security-Constrained Unit Commitment -# Copyright (C) 2020-2021, UChicago Argonne, LLC. All rights reserved. -# Released under the modified BSD license. See COPYING.md for more details. - -using Distributions - -function randomize_unit_costs!( - instance::UnitCommitmentInstance; - distribution = Uniform(0.95, 1.05), -)::Nothing - for unit in instance.units - α = rand(distribution) - unit.min_power_cost *= α - for k in unit.cost_segments - k.cost *= α - end - for s in unit.startup_categories - s.cost *= α - end - end - return -end - -function randomize_load_distribution!( - instance::UnitCommitmentInstance; - distribution = Uniform(0.90, 1.10), -)::Nothing - α = rand(distribution, length(instance.buses)) - for t in 1:instance.time - total = sum(bus.load[t] for bus in instance.buses) - den = sum( - bus.load[t] / total * α[i] for - (i, bus) in enumerate(instance.buses) - ) - for (i, bus) in enumerate(instance.buses) - bus.load[t] *= α[i] / den - end - end - return -end - -function randomize_peak_load!( - instance::UnitCommitmentInstance; - distribution = Uniform(0.925, 1.075), -)::Nothing - α = rand(distribution) - for bus in instance.buses - bus.load *= α - end - return -end - -export randomize_unit_costs!, randomize_load_distribution!, randomize_peak_load! diff --git a/src/transform/randomize/XavQiuAhm2021.jl b/src/transform/randomize/XavQiuAhm2021.jl new file mode 100644 index 0000000..adedcd2 --- /dev/null +++ b/src/transform/randomize/XavQiuAhm2021.jl @@ -0,0 +1,209 @@ +# UnitCommitment.jl: Optimization Package for Security-Constrained Unit Commitment +# Copyright (C) 2020-2021, UChicago Argonne, LLC. All rights reserved. +# Released under the modified BSD license. See COPYING.md for more details. + +""" +Methods described in: + + Xavier, Álinson S., Feng Qiu, and Shabbir Ahmed. "Learning to solve + large-scale security-constrained unit commitment problems." INFORMS + Journal on Computing 33.2 (2021): 739-756. DOI: 10.1287/ijoc.2020.0976 +""" +module XavQiuAhm2021 + +using Distributions +import ..UnitCommitmentInstance + +""" + struct Randomization + cost = Uniform(0.95, 1.05) + load_profile_mu = [...] + load_profile_sigma = [...] + load_share = Uniform(0.90, 1.10) + peak_load = Uniform(0.6 * 0.925, 0.6 * 1.075) + randomize_costs = true + randomize_load_profile = true + randomize_load_share = true + end + +Randomization method that changes: (1) production and startup costs, (2) +share of load coming from each bus, (3) peak system load, and (4) temporal +load profile, as follows: + +1. **Production and startup costs:** + For each unit `u`, the vectors `u.min_power_cost` and `u.cost_segments` + are multiplied by a constant `α[u]` sampled from the provided `cost` + distribution. If `randomize_costs` is false, skips this step. + +2. **Load share:** + For each bus `b` and time `t`, the value `b.load[t]` is multiplied by + `(β[b] * b.load[t]) / sum(β[b2] * b2.load[t] for b2 in buses)`, where + `β[b]` is sampled from the provided `load_share` distribution. If + `randomize_load_share` is false, skips this step. + +3. **Peak system load and temporal load profile:** + Sets the peak load to `ρ * C`, where `ρ` is sampled from `peak_load` and `C` + is the maximum system capacity, at any time. Also scales the loads of all + buses, so that `system_load[t+1]` becomes equal to `system_load[t] * γ[t]`, + where `γ[t]` is sampled from `Normal(load_profile_mu[t], load_profile_sigma[t])`. + + The system load for the first time period is set so that the peak load + matches `ρ * C`. If `load_profile_sigma` and `load_profile_mu` have fewer + elements than `instance.time`, wraps around. If `randomize_load_profile` + is false, skips this step. + +The default parameters were obtained based on an analysis of publicly available +bid and hourly data from PJM, corresponding to the month of January, 2017. For +more details, see Section 4.2 of the paper. +""" +Base.@kwdef struct Randomization + cost = Uniform(0.95, 1.05) + load_profile_mu::Vector{Float64} = [ + 1.0, + 0.978, + 0.98, + 1.004, + 1.02, + 1.078, + 1.132, + 1.018, + 0.999, + 1.006, + 0.999, + 0.987, + 0.975, + 0.984, + 0.995, + 1.005, + 1.045, + 1.106, + 0.981, + 0.981, + 0.978, + 0.948, + 0.928, + 0.953, + ] + load_profile_sigma::Vector{Float64} = [ + 0.0, + 0.011, + 0.015, + 0.01, + 0.012, + 0.029, + 0.055, + 0.027, + 0.026, + 0.023, + 0.013, + 0.012, + 0.014, + 0.011, + 0.008, + 0.008, + 0.02, + 0.02, + 0.016, + 0.012, + 0.014, + 0.015, + 0.017, + 0.024, + ] + load_share = Uniform(0.90, 1.10) + peak_load = Uniform(0.6 * 0.925, 0.6 * 1.075) + randomize_load_profile::Bool = true + randomize_costs::Bool = true + randomize_load_share::Bool = true +end + +function _randomize_costs( + instance::UnitCommitmentInstance, + distribution, +)::Nothing + for unit in instance.units + α = rand(distribution) + unit.min_power_cost *= α + for k in unit.cost_segments + k.cost *= α + end + for s in unit.startup_categories + s.cost *= α + end + end + return +end + +function _randomize_load_share( + instance::UnitCommitmentInstance, + distribution, +)::Nothing + α = rand(distribution, length(instance.buses)) + for t in 1:instance.time + total = sum(bus.load[t] for bus in instance.buses) + den = sum( + bus.load[t] / total * α[i] for + (i, bus) in enumerate(instance.buses) + ) + for (i, bus) in enumerate(instance.buses) + bus.load[t] *= α[i] / den + end + end + return +end + +function _randomize_load_profile( + instance::UnitCommitmentInstance, + params::Randomization, +)::Nothing + # Generate new system load + system_load = [1.0] + for t in 2:instance.time + idx = (t - 1) % length(params.load_profile_mu) + 1 + gamma = rand( + Normal(params.load_profile_mu[idx], params.load_profile_sigma[idx]), + ) + push!(system_load, system_load[t-1] * gamma) + end + capacity = sum(maximum(u.max_power) for u in instance.units) + peak_load = rand(params.peak_load) * capacity + system_load = system_load ./ maximum(system_load) .* peak_load + + # Scale bus loads to match the new system load + prev_system_load = sum(b.load for b in instance.buses) + for b in instance.buses + for t in 1:instance.time + b.load[t] *= system_load[t] / prev_system_load[t] + end + end + + return +end + +end + +""" + function randomize!( + instance::UnitCommitment.UnitCommitmentInstance, + method::XavQiuAhm2021.Randomization, + )::Nothing + +Randomize costs and loads based on the method described in XavQiuAhm2021. +""" +function randomize!( + instance::UnitCommitment.UnitCommitmentInstance, + method::XavQiuAhm2021.Randomization, +)::Nothing + if method.randomize_costs + XavQiuAhm2021._randomize_costs(instance, method.cost) + end + if method.randomize_load_share + XavQiuAhm2021._randomize_load_share(instance, method.load_share) + end + if method.randomize_load_profile + XavQiuAhm2021._randomize_load_profile(instance, method) + end + return +end + +export randomize! diff --git a/test/runtests.jl b/test/runtests.jl index 619ebd1..594d9cd 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -28,7 +28,9 @@ const ENABLE_LARGE_TESTS = ("UCJL_LARGE_TESTS" in keys(ENV)) @testset "transform" begin include("transform/initcond_test.jl") include("transform/slice_test.jl") - include("transform/randomize_test.jl") + @testset "randomize" begin + include("transform/randomize/XavQiuAhm2021_test.jl") + end end @testset "validation" begin include("validation/repair_test.jl") diff --git a/test/transform/randomize/XavQiuAhm2021_test.jl b/test/transform/randomize/XavQiuAhm2021_test.jl new file mode 100644 index 0000000..a3107c4 --- /dev/null +++ b/test/transform/randomize/XavQiuAhm2021_test.jl @@ -0,0 +1,63 @@ +# 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. + +import Random +import UnitCommitment: XavQiuAhm2021 + +using Distributions +using UnitCommitment, Cbc, JuMP + +get_instance() = UnitCommitment.read_benchmark("matpower/case118/2017-02-01") +system_load(instance) = sum(b.load for b in instance.buses) +test_approx(x, y) = @test isapprox(x, y, atol = 1e-3) + +@testset "XavQiuAhm2021" begin + @testset "cost and load share" begin + instance = get_instance() + + # Check original costs + unit = instance.units[10] + test_approx(unit.min_power_cost[1], 825.023) + test_approx(unit.cost_segments[1].cost[1], 36.659) + test_approx(unit.startup_categories[1].cost[1], 7570.42) + + # Check original load share + bus = instance.buses[1] + prev_system_load = system_load(instance) + test_approx(bus.load[1] / prev_system_load[1], 0.012) + + Random.seed!(42) + randomize!( + instance, + XavQiuAhm2021.Randomization(randomize_load_profile = false), + ) + + # Check randomized costs + test_approx(unit.min_power_cost[1], 831.977) + test_approx(unit.cost_segments[1].cost[1], 36.968) + test_approx(unit.startup_categories[1].cost[1], 7634.226) + + # Check randomized load share + curr_system_load = system_load(instance) + test_approx(bus.load[1] / curr_system_load[1], 0.013) + + # System load should not change + @test prev_system_load ≈ curr_system_load + end + + @testset "load profile" begin + instance = get_instance() + + # Check original load profile + @test round.(system_load(instance), digits = 1)[1:8] ≈ + [3059.5, 2983.2, 2937.5, 2953.9, 3073.1, 3356.4, 4068.5, 4018.8] + + Random.seed!(42) + randomize!(instance, XavQiuAhm2021.Randomization()) + + # Check randomized load profile + @test round.(system_load(instance), digits = 1)[1:8] ≈ + [4854.7, 4849.2, 4732.7, 4848.2, 4948.4, 5231.1, 5874.8, 5934.8] + end +end diff --git a/test/transform/randomize_test.jl b/test/transform/randomize_test.jl deleted file mode 100644 index a1f415e..0000000 --- a/test/transform/randomize_test.jl +++ /dev/null @@ -1,43 +0,0 @@ -# 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. - -using UnitCommitment, Cbc, JuMP - -_get_instance() = UnitCommitment.read_benchmark("matpower/case118/2017-02-01") -_total_load(instance) = sum(b.load[1] for b in instance.buses) - -@testset "randomize_unit_costs!" begin - instance = _get_instance() - unit = instance.units[10] - prev_min_power_cost = unit.min_power_cost - prev_prod_cost = unit.cost_segments[1].cost - prev_startup_cost = unit.startup_categories[1].cost - randomize_unit_costs!(instance) - @test prev_min_power_cost != unit.min_power_cost - @test prev_prod_cost != unit.cost_segments[1].cost - @test prev_startup_cost != unit.startup_categories[1].cost -end - -@testset "randomize_load_distribution!" begin - instance = _get_instance() - bus = instance.buses[1] - prev_load = instance.buses[1].load[1] - prev_total_load = _total_load(instance) - randomize_load_distribution!(instance) - curr_total_load = _total_load(instance) - @test prev_load != instance.buses[1].load[1] - @test abs(prev_total_load - curr_total_load) < 1e-3 -end - -@testset "randomize_peak_load!" begin - instance = _get_instance() - bus = instance.buses[1] - prev_total_load = _total_load(instance) - prev_share = bus.load[1] / prev_total_load - randomize_peak_load!(instance) - curr_total_load = _total_load(instance) - curr_share = bus.load[1] / prev_total_load - @test curr_total_load != prev_total_load - @test abs(curr_share - prev_share) < 1e-3 -end