From 8725973c3f49907ea46e835b984fe90289cc527c Mon Sep 17 00:00:00 2001 From: m-bossart Date: Thu, 15 Oct 2026 12:14:09 -0700 Subject: [PATCH] Fill event models from the emulator's system at simulation build Each model's build used to fill the shared event models from its own system, so the map the emulator applied at runtime depended on build order, and its attribute ids could belong to a decision model's system. The emulator's system is the single source of truth for which outages exist. Before building any model, the simulation now fills every event model from the emulator's system with POM's discover_event_attribute_device_map; each model then applies the event to the outaged devices it has. A sequence with events and no emulation model fails the build with IS.ConflictingInputsError. Outage attributes in a decision model's own system do not affect events; they may serve other formulations, such as security-constrained ones. apply_simulation_events! looks the attributes up in the emulator's system by id and checks each has the event model's type. Because the outaged devices now come from the emulator's system, an outaged device may be one no decision model has. Event-parameter state starts at its no-outage values instead of NaN, since the emulator reads its event parameters from the decision state before any event step has written them. Event-parameter keys only the emulator has get a decision state of their own, and the event-path variable writers skip an outaged device whose variable state no model writes. --- src/simulation/simulation.jl | 1 + src/simulation/simulation_events.jl | 41 +++- src/simulation/simulation_sequence.jl | 10 +- src/simulation/simulation_state.jl | 36 +++- test/test_events.jl | 262 ++++++++++++++++++++++++++ 5 files changed, 337 insertions(+), 13 deletions(-) diff --git a/src/simulation/simulation.jl b/src/simulation/simulation.jl index 70498512df..972cef96d4 100644 --- a/src/simulation/simulation.jl +++ b/src/simulation/simulation.jl @@ -600,6 +600,7 @@ function _build!( end =# + _fill_event_models!(sim) _build_decision_models!(sim) _build_emulation_model!(sim) diff --git a/src/simulation/simulation_events.jl b/src/simulation/simulation_events.jl index 80ab4d18f7..ebae8cbef2 100644 --- a/src/simulation/simulation_events.jl +++ b/src/simulation/simulation_events.jl @@ -21,11 +21,18 @@ function apply_simulation_events!(simulation::Simulation) # updated every step regardless, or a running outage would freeze the moment the # condition stopped holding and never recover. may_start = _condition_holds(simulation, event_model) + event_type = POM.get_event_type(event_model) for (attribute_id, device_type_maps) in POM.get_attribute_device_map(event_model) + attribute = PSY.get_supplemental_attribute(sys, attribute_id) + attribute isa event_type || error( + "Supplemental attribute $attribute_id of the emulator's system is a \ + $(typeof(attribute)), but the $event_type event model was filled with \ + it. Event models must be filled from the emulator's system.", + ) apply_event_step!( simulation, event_model, - PSY.get_supplemental_attribute(sys, attribute_id), + attribute, device_type_maps; may_start = may_start, ) @@ -34,6 +41,38 @@ function apply_simulation_events!(simulation::Simulation) return end +""" +Fill every event model of the simulation from the emulator's system, before any model is +built. Events take their outages from the emulator's system only; each model then applies +the event to the outaged devices it has. Outage attributes in a decision model's own system +do not affect events; they may serve other formulations, such as security-constrained ones. + +Throws `IS.ConflictingInputsError` when the simulation has event models but no emulation +model. +""" +function _fill_event_models!(simulation::Simulation) + events = get_events(get_sequence(simulation)) + isempty(events) && return + models = get_models(simulation) + em_model = get_emulation_model(models) + if isnothing(em_model) + throw( + IS.ConflictingInputsError( + "The simulation sequence has event models but the simulation has no \ + emulation model. Outage events are drawn and applied by the emulation \ + model from its system: add an emulation model or remove the events from \ + the SimulationSequence.", + ), + ) + end + em_sys = get_system(em_model) + for event_model in events + outages = POM.discover_event_attribute_device_map(event_model, em_sys) + copy!(POM.get_attribute_device_map(event_model), outages) + end + return +end + """ Resolve one input an event condition declared. POM says what it needs; this says where it comes from. diff --git a/src/simulation/simulation_sequence.jl b/src/simulation/simulation_sequence.jl index 3f6d3cea5b..a60d79fbb6 100644 --- a/src/simulation/simulation_sequence.jl +++ b/src/simulation/simulation_sequence.jl @@ -201,11 +201,11 @@ end """ Attach every event model to the template of every model in the simulation. -Discovery and validation happen in PowerOperationsModels at build time: attaching an -event model to a template is what makes `build!` find the system's outage attributes, -populate the event model's `attribute_device_map`, and distribute the event to the -matching `DeviceModel`s. This function used to do that work itself, per model, against -a map keyed by model name. +Every model shares the same event model objects. At simulation build, +`_fill_event_models!` fills each event model's `attribute_device_map` from the emulator's +system; each model's build then distributes the event to its matching `DeviceModel`s and +applies it to the outaged devices that model has. Outage attributes in a decision model's +own system do not affect events. """ function _attach_events!( models::SimulationModels, diff --git a/src/simulation/simulation_state.jl b/src/simulation/simulation_state.jl index 879808b173..22b0885640 100644 --- a/src/simulation/simulation_state.jl +++ b/src/simulation/simulation_state.jl @@ -79,7 +79,11 @@ function _get_state_params(models::SimulationModels, simulation_step::Dates.Mill field_containers = getfield(container, type) for key in keys(field_containers) !should_write_resulting_value(key) && continue - if !haskey(params, key) + if !haskey(params, key) && _is_event_parameter_key(key) + # The emulator reads its event parameters from the decision state, so + # a device type only the emulator applies events to still needs it. + params[key] = (horizon = simulation_step, resolution = model_resolution) + elseif !haskey(params, key) @warn "New parameter $key found in emulator only; \ not registered by any decision model" else @@ -95,6 +99,9 @@ function _get_state_params(models::SimulationModels, simulation_step::Dates.Mill return params end +_is_event_parameter_key(::OptimizationContainerKey) = false +_is_event_parameter_key(::ParameterKey{<:EventParameter}) = true + _state_models(models::SimulationModels) = _state_models(models, get_emulation_model(models)) _state_models(models::SimulationModels, ::Nothing) = get_decision_models(models) _state_models(models::SimulationModels, em::EmulationModel) = @@ -161,6 +168,17 @@ function _get_state_columns(models::SimulationModels) return columns end +# Event parameters are exogenous state the emulator reads before any event step writes it. +# Starting them at their no-outage values covers the components no decision model writes, +# such as a device only the emulator has; every other key starts unset. +_initial_state_value(::OptimizationContainerKey) = NaN +_initial_state_value(::ParameterKey{<:EventParameter}) = 0.0 +_initial_state_value(::ParameterKey{AvailableStatusParameter}) = 1.0 + +# The outaged devices of an event are those of the emulator's system, but a variable's state +# holds only the components of the models that write it, which need not include all of them. +_has_component(data::InMemoryDataset, name::String) = name in axes(data.values)[1] + function _initialize_decision_states!( sim_state::SimulationState, simulation_initial_time::Dates.DateTime, @@ -174,7 +192,7 @@ function _initialize_decision_states!( states, key, InMemoryDataset( - NaN, + _initial_state_value(key), simulation_initial_time, key_params.resolution, Int(simulation_step / key_params.resolution), @@ -230,11 +248,9 @@ function _set_system_states!( system_states = get_system_states(sim_state) min_res = minimum([v.resolution for v in values(params)]) for key in state_keys - set_dataset!( - system_states, - key, - make_system_state(simulation_initial_time, min_res, columns[key]), - ) + dataset = make_system_state(simulation_initial_time, min_res, columns[key]) + fill!(dataset.values.data, _initial_state_value(key)) + set_dataset!(system_states, key, dataset) end return end @@ -541,6 +557,7 @@ function update_decision_state!( event_occurence_index = find_aligned_timestamp_index(event_occurrence_data, simulation_time) for name in column_names + _has_component(state_data, name) || continue if event_occurrence_data.values[name, event_occurence_index] > 0.0 state_data.values[name, (state_data_index + 1):end] .= MISSING_INITIAL_CONDITIONS_TIME_COUNT @@ -566,6 +583,7 @@ function update_decision_state!( event_occurence_index = find_aligned_timestamp_index(event_occurrence_data, simulation_time) for name in column_names + _has_component(state_data, name) || continue if event_occurrence_data.values[name, event_occurence_index] > 0.0 for (time_off, ix) in enumerate((state_data_index + 1):length(state_data.values[name, :])) @@ -593,6 +611,7 @@ function update_decision_state!( event_occurence_index = find_aligned_timestamp_index(event_occurrence_data, simulation_time) for name in column_names + _has_component(state_data, name) || continue if event_occurrence_data.values[name, event_occurence_index] > 0.0 state_data.values[name, (state_data_index + 1):end] .= 0.0 end @@ -1073,6 +1092,7 @@ function update_system_state!( current_status_values = get_last_recorded_value(current_status_data) set_update_timestamp!(system_dataset, simulation_time) for name in column_names + _has_component(current_status_data, name) || continue if event_occurrence_values[name] > 0.0 current_status_data.values[name, 1] = 0.0 end @@ -1099,6 +1119,7 @@ function update_system_state!( current_status_values = get_last_recorded_value(current_status_data) set_update_timestamp!(system_dataset, simulation_time) for name in column_names + _has_component(current_status_data, name) || continue if event_occurrence_values[name] > 0.0 current_status_data.values[name, 1] = MISSING_INITIAL_CONDITIONS_TIME_COUNT end @@ -1125,6 +1146,7 @@ function update_system_state!( current_status_values = get_last_recorded_value(current_status_data) set_update_timestamp!(system_dataset, simulation_time) for name in column_names + _has_component(current_status_data, name) || continue if event_occurrence_values[name] > 0.0 current_status_data.values[name, 1] = 0.0 end diff --git a/test/test_events.jl b/test/test_events.jl index dc3c662a04..3cfb08a3b5 100644 --- a/test/test_events.jl +++ b/test/test_events.jl @@ -550,3 +550,265 @@ end test_reactive_power = true, ) end + +# c_sys5_events without the thermal units in `removed_gens`, and with a +# GeometricDistributionForcedOutage, drawn with certainty, on each of `outage_gens`. +function events_system(; + outage_gens = String[], + removed_gens = String[], + horizon_interval = nothing, +) + sys = build_system(PSITestSystems, "c_sys5_events") + for name in removed_gens + remove_component!(sys, get_component(ThermalStandard, sys, name)) + end + isnothing(horizon_interval) || transform_single_time_series!(sys, horizon_interval...) + for name in outage_gens + add_supplemental_attribute!( + sys, + get_component(ThermalStandard, sys, name), + PSY.GeometricDistributionForcedOutage(; + mean_time_to_recovery = 3.0, + outage_transition_probability = 1.0, + ), + ) + end + return sys +end + +function events_simulation( + sys_d1, + sys_d2, + sys_em, + event_model; + rebuild_model = false, +) + templates = map(1:2) do _ + template = get_template_basic_uc_simulation() + set_network_model!(template, NetworkModel(CopperPlateNetworkModel)) + set_device_model!(template, Line, StaticBranch) + set_service_model!(template, ServiceModel(OnlineReserve{ReserveUp}, RangeReserve)) + template + end + decision_models = [ + DecisionModel( + templates[1], + sys_d1; + name = "D1", + initialize_model = false, + optimizer = HiGHS_optimizer_small_gap, + rebuild_model = rebuild_model, + ), + DecisionModel( + templates[2], + sys_d2; + name = "D2", + initialize_model = false, + optimizer = HiGHS_optimizer_small_gap, + rebuild_model = rebuild_model, + ), + ] + if isnothing(sys_em) + models = SimulationModels(; decision_models = decision_models) + else + template_em = get_template_nomin_ed_simulation(CopperPlateNetworkModel) + set_device_model!(template_em, Line, StaticBranch) + models = SimulationModels(; + decision_models = decision_models, + emulation_model = EmulationModel( + template_em, + sys_em; + name = "EM", + optimizer = HiGHS_optimizer_small_gap, + rebuild_model = rebuild_model, + ), + ) + end + sequence = SimulationSequence(; + models = models, + ini_cond_chronology = InterProblemChronology(), + events = [event_model], + ) + return Simulation(; + name = "events", + steps = 1, + models = models, + sequence = sequence, + simulation_folder = mktempdir(; cleanup = true), + ) +end + +preset_outage_event(outage_time) = + EventModel(GeometricDistributionForcedOutage, POM.PresetTimeCondition([outage_time])) + +function event_device_names(model) + status = IOM.get_parameter_array( + IOM.get_optimization_container(model), + AvailableStatusParameter(), + ThermalStandard, + ) + return Set(axes(status)[1]) +end + +@testset "Outages are applied from the emulator's system when attribute ids diverge" begin + outage_time = DateTime("2024-01-01T18:00:00") + sys_em = build_system(PSITestSystems, "c_sys5_events") + bus = first(get_components(ACBus, sys_em)) + add_supplemental_attribute!( + sys_em, + bus, + PSY.GeographicInfo(; + geo_json = Dict("type" => "Point", "coordinates" => [0.0, 0.0]), + ), + ) + res = run_events_simulation(; + sys_emulator = sys_em, + networks = repeat([CopperPlateNetworkModel], 3), + optimizers = repeat([HiGHS_optimizer_small_gap], 3), + outage_time = outage_time, + outage_length = 3.0, + uc_formulation = "basic", + ed_formulation = "nomin", + feedforward = false, + in_memory = true, + ) + test_event_results(; + res = res, + outage_time = outage_time, + outage_length = 3.0, + expected_power_recovery = DateTime("2024-01-01T22:00:00"), + expected_on_variable_recovery = DateTime("2024-01-01T22:00:00"), + ) +end + +@testset "Every decision model sees the emulator's outages, and only those" begin + # The emulator's system alone decides which outages exist; both of its outages are + # drawn at `outage_time`. D1's own outage on Brighton does not affect events. + outage_time = DateTime("2024-01-01T18:00:00") + outaged = ["Alta", "Park City"] + sys_d1 = events_system(; + outage_gens = ["Brighton"], + horizon_interval = (Day(2), Day(1)), + ) + sys_d2 = events_system(; horizon_interval = (Hour(4), Hour(1))) + sys_em = events_system(; outage_gens = outaged) + event_model = preset_outage_event(outage_time) + sim = events_simulation(sys_d1, sys_d2, sys_em, event_model) + @test build!(sim; console_level = Logging.Error) == SimulationBuildStatus.BUILT + + @test POM.get_attribute_device_map(event_model) == Dict( + IS.get_id(outage) => + Dict( + ThermalStandard => + Set(get_name.(get_associated_components(sys_em, outage))), + ) + for outage in get_supplemental_attributes(GeometricDistributionForcedOutage, sys_em) + ) + for model in PSI.get_decision_models(PSI.get_models(sim)) + @test event_device_names(model) == Set(outaged) + end + + @test execute!(sim; in_memory = true) == IOM.RunStatus.SUCCESSFULLY_FINALIZED + d2 = get_decision_problem_results(SimulationResults(sim; ignore_status = true), "D2") + read_d2(name) = read_realized_variable(d2, name; table_format = TableFormat.WIDE) + status = read_d2("AvailableStatusParameter__ThermalStandard") + p = read_d2("ActivePowerVariable__ThermalStandard") + @test Set(names(status, Not(:DateTime))) == Set(outaged) + out = outage_time + Hour(1) .<= status.DateTime .<= outage_time + Hour(3) + for name in outaged + @test all(status[status.DateTime .<= outage_time, name] .== 1.0) + @test all(status[out, name] .== 0.0) + @test isapprox(p[out, name], zeros(count(out)); atol = 1e-5) + end +end + +# Run a simulation whose emulator outages `outaged` while the decision systems lack +# `emulator_only`, and return the realized results of the emulator and of D2. +function run_emulator_only_outage(outaged, emulator_only, outage_time) + sys_d1 = events_system(; + removed_gens = emulator_only, + horizon_interval = (Day(2), Day(1)), + ) + sys_d2 = events_system(; + removed_gens = emulator_only, + horizon_interval = (Hour(4), Hour(1)), + ) + sys_em = events_system(; outage_gens = outaged) + sim = events_simulation(sys_d1, sys_d2, sys_em, preset_outage_event(outage_time)) + @test build!(sim; console_level = Logging.Error) == SimulationBuildStatus.BUILT + @test execute!(sim; in_memory = true) == IOM.RunStatus.SUCCESSFULLY_FINALIZED + res = SimulationResults(sim; ignore_status = true) + return get_emulation_problem_results(res), get_decision_problem_results(res, "D2") +end + +read_wide(res, name) = read_realized_variable(res, name; table_format = TableFormat.WIDE) + +# The device is out for the three hours after `outage_time` and available before. +function test_outaged(res, name, outage_time) + status = read_wide(res, "AvailableStatusParameter__ThermalStandard") + p = read_wide(res, "ActivePowerVariable__ThermalStandard") + out = outage_time + Hour(1) .<= status.DateTime .<= outage_time + Hour(3) + @test all(status[status.DateTime .<= outage_time, name] .== 1.0) + @test all(status[out, name] .== 0.0) + @test isapprox(p[out, name], zeros(count(out)); atol = 1e-5) + return +end + +@testset "An outage on a device only the emulator has" begin + outage_time = DateTime("2024-01-01T18:00:00") + em, d2 = run_emulator_only_outage(["Alta", "Park City"], ["Park City"], outage_time) + test_outaged(em, "Park City", outage_time) + test_outaged(em, "Alta", outage_time) + test_outaged(d2, "Alta", outage_time) + status = read_wide(d2, "AvailableStatusParameter__ThermalStandard") + @test names(status, Not(:DateTime)) == ["Alta"] + @test !("Park City" in names(read_wide(d2, "ActivePowerVariable__ThermalStandard"))) +end + +@testset "An outage of a device type no decision model has event parameters for" begin + outage_time = DateTime("2024-01-01T18:00:00") + em, d2 = run_emulator_only_outage(["Park City"], ["Park City"], outage_time) + test_outaged(em, "Park City", outage_time) + @test !( + "AvailableStatusParameter__ThermalStandard" in + list_parameter_names(d2) + ) +end + +@testset "Events without an emulation model fail the build" begin + sys_d1 = events_system(; outage_gens = ["Alta"], horizon_interval = (Day(2), Day(1))) + sys_d2 = events_system(; outage_gens = ["Alta"], horizon_interval = (Hour(4), Hour(1))) + event_model = preset_outage_event(DateTime("2024-01-01T18:00:00")) + sim = events_simulation(sys_d1, sys_d2, nothing, event_model) + @test_throws IS.ConflictingInputsError build!(sim; console_level = Logging.Error) + @test PSI.get_simulation_build_status(sim) == SimulationBuildStatus.FAILED + @test isempty(POM.get_attribute_device_map(event_model)) +end + +@testset "Rebuilt models with events (rebuild_model = true)" begin + outage_time = DateTime("2024-01-01T18:00:00") + sys_d1 = events_system(; horizon_interval = (Day(2), Day(1))) + sys_d2 = events_system(; horizon_interval = (Hour(4), Hour(1))) + sys_em = events_system(; outage_gens = ["Alta"]) + event_model = preset_outage_event(outage_time) + sim = events_simulation( + sys_d1, + sys_d2, + sys_em, + event_model; + rebuild_model = true, + ) + @test all( + IOM.get_rebuild_model(IOM.get_settings(model)) for + model in PSI.get_decision_models(PSI.get_models(sim)) + ) + @test build!(sim; console_level = Logging.Error) == SimulationBuildStatus.BUILT + @test execute!(sim; in_memory = true) == IOM.RunStatus.SUCCESSFULLY_FINALIZED + test_event_results(; + res = SimulationResults(sim; ignore_status = true), + outage_time = outage_time, + outage_length = 3.0, + expected_power_recovery = DateTime("2024-01-01T22:00:00"), + expected_on_variable_recovery = DateTime("2024-01-01T22:00:00"), + ) +end