Originally Contributed by: Clayton Barrows
PowerSimulations.jl supports simulations that consist of sequential optimization problems where results from previous problems inform subsequent problems in a variety of ways. This example demonstrates some of these capabilities to represent electricity market clearing.
First, let's create Systems to represent the Day-Ahead and Real-Time market clearing
process with hourly, and 5-minute time series data, respectively.
using PowerSystems
using PowerSimulations
using PowerSystemCaseBuilder
using Dates
using DataFrames
using Cbc # mip solver
solver = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 1, "ratioGap" => 0.5)
using Ipopt # solver that supports duals
ipopt_solver = optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)
MathOptInterface.OptimizerWithAttributes(Ipopt.Optimizer, Pair{MathOptInterface.AbstractOptimizerAttribute,Any}[MathOptInterface.RawParameter("print_level") => 0])
The five bus system data here includes hourly day-ahead data, 5-minute real-time market data, and 6-second actual data. We'll only use the hourly and 5-minute data for the example simulations below, but the 6-second data is included for future development.
sys_DA = build_system(SIIPExampleSystems, "5_bus_matpower_DA")
sys_RT = build_system(SIIPExampleSystems, "5_bus_matpower_RT")
┌ Info: Building new system 5_bus_matpower_DA from raw data └ sys_descriptor.raw_data = "/Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/matpower/case5_re_uc.m" [ Info: extending matpower format with data: areas 1x3 [ Info: extending matpower format with data: gen_name 7x4 [ Info: extending matpower format by appending matrix "gen_name" in to "gen" [ Info: reversing the orientation of branch 6 (4, 3) to be consistent with other parallel branches [ Info: the voltage setpoint on generator 4 does not match the value at bus 4 [ Info: the voltage setpoint on generator 1 does not match the value at bus 1 [ Info: the voltage setpoint on generator 5 does not match the value at bus 10 [ Info: the voltage setpoint on generator 2 does not match the value at bus 1 [ Info: the voltage setpoint on generator 3 does not match the value at bus 3 [ Info: removing 1 cost terms from generator 4: [4000.0, 0.0] [ Info: removing 1 cost terms from generator 1: [1400.0, 0.0] [ Info: removing 1 cost terms from generator 5: [1000.0, 0.0] [ Info: removing 1 cost terms from generator 2: [1500.0, 0.0] [ Info: removing 3 cost terms from generator 6: Float64[] [ Info: removing 3 cost terms from generator 7: Float64[] [ Info: removing 1 cost terms from generator 3: [3000.0, 0.0] ┌ Info: Constructing System from Power Models │ data["name"] = "nesta_case5_pjm" └ data["source_type"] = "matpower" [ Info: Reading bus data [ Info: Reading generator data [ Info: Reading branch data ┌ Warning: Rate 1200.0 MW for bus2-bus3-i_4 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus1-bus2-i_1 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus1-bus4-i_2 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus4-bus5-i_7 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 1000.0 MW for bus1-bus5-i_3 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 [ Info: Reading branch data [ Info: Reading DC Line data [ Info: Reading storage data ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch2 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch2 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 [ Info: Serialized time series data to /Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/serialized_system/ForecastOnly/5_bus_matpower_DA_time_series_storage.h5. [ Info: Serialized System to /Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/serialized_system/ForecastOnly/5_bus_matpower_DA.json ┌ Info: Building new system 5_bus_matpower_RT from raw data └ sys_descriptor.raw_data = "/Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/matpower/case5_re_uc.m" [ Info: extending matpower format with data: areas 1x3 [ Info: extending matpower format with data: gen_name 7x4 [ Info: extending matpower format by appending matrix "gen_name" in to "gen" [ Info: reversing the orientation of branch 6 (4, 3) to be consistent with other parallel branches [ Info: the voltage setpoint on generator 4 does not match the value at bus 4 [ Info: the voltage setpoint on generator 1 does not match the value at bus 1 [ Info: the voltage setpoint on generator 5 does not match the value at bus 10 [ Info: the voltage setpoint on generator 2 does not match the value at bus 1 [ Info: the voltage setpoint on generator 3 does not match the value at bus 3 [ Info: removing 1 cost terms from generator 4: [4000.0, 0.0] [ Info: removing 1 cost terms from generator 1: [1400.0, 0.0] [ Info: removing 1 cost terms from generator 5: [1000.0, 0.0] [ Info: removing 1 cost terms from generator 2: [1500.0, 0.0] [ Info: removing 3 cost terms from generator 6: Float64[] [ Info: removing 3 cost terms from generator 7: Float64[] [ Info: removing 1 cost terms from generator 3: [3000.0, 0.0] ┌ Info: Constructing System from Power Models │ data["name"] = "nesta_case5_pjm" └ data["source_type"] = "matpower" [ Info: Reading bus data [ Info: Reading generator data [ Info: Reading branch data ┌ Warning: Rate 1200.0 MW for bus2-bus3-i_4 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus1-bus2-i_1 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus1-bus4-i_2 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 200.0 MW for bus4-bus5-i_7 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 ┌ Warning: Rate 1000.0 MW for bus1-bus5-i_3 is larger than the max expected in the range of (min = 134.0, max = 145.0). └ @ PowerSystems ~/.julia/packages/PowerSystems/N2l8o/src/utils/IO/branchdata_checks.jl:148 [ Info: Reading branch data [ Info: Reading DC Line data [ Info: Reading storage data ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch2 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch1 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 ┌ Warning: no component category=Generator name=HydroDispatch3 └ @ InfrastructureSystems ~/.julia/packages/InfrastructureSystems/lYELp/src/system_data.jl:448 [ Info: Serialized time series data to /Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/serialized_system/ForecastOnly/5_bus_matpower_RT_time_series_storage.h5. [ Info: Serialized System to /Users/cbarrows/.julia/packages/PowerSystemCaseBuilder/1xViZ/data/serialized_system/ForecastOnly/5_bus_matpower_RT.json
Base Power: 100.0
Num components: 30
| ConcreteType | SuperTypes | Count | |
|---|---|---|---|
| String | String | Int64 | |
| 1 | Arc | Topology <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 6 |
| 2 | Area | AggregationTopology <: Topology <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 1 |
| 3 | Bus | Topology <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 5 |
| 4 | Line | ACBranch <: Branch <: Device <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 5 |
| 5 | LoadZone | AggregationTopology <: Topology <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 1 |
| 6 | PhaseShiftingTransformer | ACBranch <: Branch <: Device <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 2 |
| 7 | PowerLoad | StaticLoad <: ElectricLoad <: StaticInjection <: Device <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 3 |
| 8 | RenewableDispatch | RenewableGen <: Generator <: StaticInjection <: Device <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 2 |
| 9 | ThermalStandard | ThermalGen <: Generator <: StaticInjection <: Device <: Component <: InfrastructureSystemsComponent <: InfrastructureSystemsType <: Any | 5 |
Components with time series data: 5
Total StaticTimeSeries: 5
Total Forecasts: 5
Resolution: 60 minutes
First initial time: 2020-01-01T00:00:00
Last initial time: 2020-03-24T12:00:00
Horizon: 12
Interval: 60 minutes
Forecast window count: 2005
OperationsProblemTemplates¶template_uc = template_unit_commitment()
template_ed = template_economic_dispatch()
Operations Problem Specification
============================================
Transmission: CopperPlatePowerModel
============================================
Devices Models:
Type: ThermalStandard
Formulation: ThermalDispatch
Type: HydroDispatch
Formulation: HydroDispatchRunOfRiver
Type: PowerLoad
Formulation: StaticPowerLoad
Type: RenewableFix
Formulation: FixedOutput
Type: RenewableDispatch
Formulation: RenewableFullDispatch
Type: HydroEnergyReservoir
Formulation: HydroDispatchRunOfRiver
Type: InterruptibleLoad
Formulation: InterruptiblePowerLoad
============================================
Branches Models:
Type: Line
Formulation: StaticBranch
Type: TapTransformer
Formulation: StaticBranch
Type: Transformer2W
Formulation: StaticBranch
Type: HVDCLine
Formulation: HVDCDispatch
============================================
Services Models:
Type: VariableReserve{ReserveDown}
Formulation: RangeReserve
Type: VariableReserve{ReserveUp}
Formulation: RangeReserve
============================================
problems = SimulationProblems(
UC = OperationsProblem(
template_uc,
sys_DA,
optimizer = solver,
balance_slack_variables = true,
),
ED = OperationsProblem(
template_ed,
sys_RT,
optimizer = ipopt_solver,
constraint_duals = [:CopperPlateBalance],
),
)
feedforward_chronologies = Dict(("UC" => "ED") => Synchronize(periods = 24))
feedforward = Dict(
("ED", :devices, :ThermalStandard) => SemiContinuousFF(
binary_source_problem = ON,
affected_variables = [ACTIVE_POWER],
),
)
#cache = Dict("UC" => [TimeStatusChange(ThermalStandard, PSI.ON)])
intervals = Dict("UC" => (Hour(24), Consecutive()), "ED" => (Hour(1), Consecutive()))
DA_RT_sequence = SimulationSequence(
problems = problems,
intervals = intervals,
ini_cond_chronology = InterProblemChronology(),
feedforward_chronologies = feedforward_chronologies,
feedforward = feedforward,
)
Feed Forward Chronology
-----------------------
ED: SemiContinuousFF -> ThermalStandard
UC--┐ from : On
|
┌----┬----┬----┬----┬----┼----┬----┬----┬----┬----┬----┐
| | | | | | | | | | | |
| | | | | | | | | | | |
└─ED └─ED └─ED └─ED └─ED └─ED └─ED └─ED └─ED └─ED └─ED └─ED ... (x24) to : ["P"]
Initial Condition Chronology
----------------------------
1
|
|
2 --> 2 ... (x24)
Simulation¶file_path = mkpath(joinpath(".", "5-bus-simulation"))
sim = Simulation(
name = "5bus-test",
steps = 1,
problems = problems,
sequence = DA_RT_sequence,
simulation_folder = file_path,
)
Simulation()
build!(sim)
┌ Info: │ ────────────────────────────────────────────────────────────────────────────── │ Time Allocations │ ────────────────────── ─────────────────────── │ Tot / % measured: 1.78s / 100% 148MiB / 100% │ │ Section ncalls time %tot avg alloc %tot avg │ ────────────────────────────────────────────────────────────────────────────── │ Build Simulation 1 1.78s 100% 1.78s 148MiB 100% 148MiB │ Build Problems 1 1.77s 100% 1.77s 148MiB 100% 148MiB │ Problem UC 1 1.21s 68.3% 1.21s 117MiB 79.2% 117MiB │ ThermalStandard 1 521ms 29.3% 521ms 41.0MiB 27.7% 41.0MiB │ Services 1 96.8ms 5.45% 96.8ms 10.5MiB 7.06% 10.5MiB │ Build pre-step 1 29.4ms 1.65% 29.4ms 2.79MiB 1.89% 2.79MiB │ RenewableDis... 1 4.12ms 0.23% 4.12ms 962KiB 0.63% 962KiB │ PowerLoad 1 2.98ms 0.17% 2.98ms 321KiB 0.21% 321KiB │ CopperPlateP... 1 747μs 0.04% 747μs 454KiB 0.30% 454KiB │ Objective 1 78.3μs 0.00% 78.3μs 109KiB 0.07% 109KiB │ HydroDispatch 1 59.2μs 0.00% 59.2μs 4.58KiB 0.00% 4.58KiB │ HydroEnergyR... 1 44.0μs 0.00% 44.0μs 4.70KiB 0.00% 4.70KiB │ RenewableFix 1 40.8μs 0.00% 40.8μs 4.58KiB 0.00% 4.58KiB │ Interruptibl... 1 17.7μs 0.00% 17.7μs 4.63KiB 0.00% 4.63KiB │ Line 1 2.40μs 0.00% 2.40μs 64.0B 0.00% 64.0B │ HVDCLine 1 1.44μs 0.00% 1.44μs 64.0B 0.00% 64.0B │ Transformer2W 1 1.22μs 0.00% 1.22μs 64.0B 0.00% 64.0B │ TapTransformer 1 1.11μs 0.00% 1.11μs 64.0B 0.00% 64.0B │ Problem ED 1 560ms 31.5% 560ms 30.6MiB 20.6% 30.6MiB │ PowerLoad 1 3.56ms 0.20% 3.56ms 113KiB 0.07% 113KiB │ RenewableDis... 1 1.85ms 0.10% 1.85ms 250KiB 0.17% 250KiB │ ThermalStandard 1 1.78ms 0.10% 1.78ms 915KiB 0.60% 915KiB │ Build pre-step 1 416μs 0.02% 416μs 61.7KiB 0.04% 61.7KiB │ CopperPlateP... 1 117μs 0.01% 117μs 54.3KiB 0.04% 54.3KiB │ Services 1 98.3μs 0.01% 98.3μs 10.4KiB 0.01% 10.4KiB │ RenewableFix 1 47.1μs 0.00% 47.1μs 4.58KiB 0.00% 4.58KiB │ HydroEnergyR... 1 33.2μs 0.00% 33.2μs 4.70KiB 0.00% 4.70KiB │ HydroDispatch 1 26.7μs 0.00% 26.7μs 4.58KiB 0.00% 4.58KiB │ Interruptibl... 1 20.4μs 0.00% 20.4μs 4.63KiB 0.00% 4.63KiB │ Objective 1 13.4μs 0.00% 13.4μs 5.22KiB 0.00% 5.22KiB │ Line 1 2.86μs 0.00% 2.86μs 64.0B 0.00% 64.0B │ HVDCLine 1 1.80μs 0.00% 1.80μs 64.0B 0.00% 64.0B │ TapTransformer 1 1.34μs 0.00% 1.34μs 64.0B 0.00% 64.0B │ Transformer2W 1 1.31μs 0.00% 1.31μs 64.0B 0.00% 64.0B │ Initialize Simul... 1 1.50ms 0.08% 1.50ms 19.7KiB 0.01% 19.7KiB │ Serializing Simu... 1 757μs 0.04% 757μs 16.0KiB 0.01% 16.0KiB │ Assign FeedForward 1 35.4μs 0.00% 35.4μs 3.16KiB 0.00% 3.16KiB │ Check Steps 1 25.0μs 0.00% 25.0μs 4.45KiB 0.00% 4.45KiB └ ──────────────────────────────────────────────────────────────────────────────
BuildStatus.BUILT = 0
execute!(sim, enable_progress_bar = false)
Presolve 701 (-2599) rows, 1589 (-5971) columns and 2930 (-11098) elements
Perturbing problem by 0.001% of 10 - largest nonzero change 0.00010943053 ( 54.949531%) - largest zero change 0.00010938892
0 Obj 0.0082926214 Primal inf 370.67927 (425)
89 Obj 13060.303 Primal inf 182.89395 (301)
178 Obj 24299.651 Primal inf 100.37047 (232)
267 Obj 30950.67 Primal inf 56.412648 (179)
356 Obj 31158.922 Primal inf 23.524935 (120)
445 Obj 31571.851 Primal inf 5.3320629 (47)
498 Obj 31907.354
Optimal - objective value 31907.325
After Postsolve, objective 31907.325, infeasibilities - dual 42318.909 (240), primal 0 (0)
Presolved model was optimal, full model needs cleaning up
0 Obj 31907.325
Optimal - objective value 31907.325
Optimal objective 31907.32482 - 498 iterations time 0.012, Presolve 0.01
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit https://github.com/coin-or/Ipopt
******************************************************************************
RunStatus.SUCCESSFUL = 0
# Results
First we can load the result metadata
results = SimulationResults(sim);
uc_results = get_problem_results(results, "UC")
ed_results = get_problem_results(results, "ED");
Then we can read and examine the results of interest
prices = read_dual(ed_results, :CopperPlateBalance)
or if we want to look at the realized values
read_realized_duals(ed_results)[:CopperPlateBalance]
note that in this simulation the prices are all equal to the balance slack penalty value of $100000/MWh because there is unserved energy in the result
This notebook was generated using Literate.jl.