diff --git a/Project.toml b/Project.toml index fa73c032..25922a2e 100644 --- a/Project.toml +++ b/Project.toml @@ -15,6 +15,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" PowerNetworkMatrices = "bed98974-b02a-5e2f-9fe0-a103f5c450dd" PowerSystems = "bcd98974-b02a-5e2f-9ee0-a103f5c450dd" +PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d" ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" Serialization = "9e88b42a-f829-5b0c-bbe9-9e923198166b" @@ -42,6 +43,7 @@ JuMP = "^1.28" PowerSystems = "^5.10" PowerFlows = "0.21, 0.22" PowerNetworkMatrices = "^0.24" +PrecompileTools = "1" PrettyTables = "3" ProgressMeter = "1.11.0" TimerOutputs = "~0.5" diff --git a/src/PowerOperationsModels.jl b/src/PowerOperationsModels.jl index 7e0d2551..35054d70 100644 --- a/src/PowerOperationsModels.jl +++ b/src/PowerOperationsModels.jl @@ -14,6 +14,7 @@ import PowerNetworkMatrices import PowerNetworkMatrices: PTDF, VirtualPTDF import PowerSystems import PowerSystems: get_component +import PrecompileTools import PrettyTables import ProgressMeter import Serialization @@ -973,4 +974,8 @@ export PostContingencyBranchRatingTimeSeriesParameter # Network model capability traits (IOM default + POM overrides in core/network_formulations.jl) export branches_modeled +# Must be the last include: it runs build!/solve! at precompile time and +# needs the post-include IOM imports (set_status! etc.) already evaluated. +include("precompile_workload.jl") + end diff --git a/src/precompile_workload.jl b/src/precompile_workload.jl new file mode 100644 index 00000000..ad2a05ff --- /dev/null +++ b/src/precompile_workload.jl @@ -0,0 +1,318 @@ +# Executes build!/solve! at precompile time, so this file must stay the last +# include of the module: the IOM method imports after the main include block +# (set_status! etc.) have to be evaluated first. + +function _build_precompile_system() + sys = PSY.System(100.0; time_series_in_memory = true) + # One area per bus so the AreaPTDF template has a non-degenerate area graph. + area1 = PSY.Area(; name = "area1") + area2 = PSY.Area(; name = "area2") + PSY.add_component!(sys, area1) + PSY.add_component!(sys, area2) + bus1 = PSY.ACBus(; + number = 1, + name = "bus1", + available = true, + bustype = PSY.ACBusTypes.REF, + angle = 0.0, + magnitude = 1.0, + voltage_limits = (0.9, 1.1), + base_voltage = 230.0, + area = area1, + ) + bus2 = PSY.ACBus(; + number = 2, + name = "bus2", + available = true, + bustype = PSY.ACBusTypes.PV, + angle = 0.0, + magnitude = 1.0, + voltage_limits = (0.9, 1.1), + base_voltage = 230.0, + area = area2, + ) + PSY.add_component!(sys, bus1) + PSY.add_component!(sys, bus2) + arc = PSY.Arc(; from = bus1, to = bus2) + PSY.add_component!(sys, arc) + line = PSY.Line(; + name = "line1", + available = true, + active_power_flow = 0.0, + reactive_power_flow = 0.0, + arc = arc, + r = 0.01, + x = 0.05, + b = (from = 0.0, to = 0.0), + rating = 5.0, + angle_limits = (min = -1.57, max = 1.57), + ) + PSY.add_component!(sys, line) + thermal1 = PSY.ThermalMultiStart(; + name = "thermal1", + available = true, + status = true, + bus = bus1, + active_power = 0.5, + reactive_power = 0.0, + rating = 1.0, + prime_mover_type = PSY.PrimeMovers.CC, + fuel = PSY.ThermalFuels.NATURAL_GAS, + active_power_limits = (min = 0.2, max = 1.0), + reactive_power_limits = nothing, + # Must stay below |max - min| / 60 (pu/min at hourly resolution) or + # _get_ramp_constraint_devices drops the device and no RampConstraint + # is compiled. Same logic for thermal2. + ramp_limits = (up = 0.005, down = 0.005), + power_trajectory = (startup = 0.25, shutdown = 0.25), + time_limits = (up = 2.0, down = 2.0), + start_time_limits = (hot = 2.0, warm = 4.0, cold = 12.0), + start_types = 3, + operation_cost = PSY.ThermalGenerationCost(; + variable = PSY.CostCurve(PSY.LinearCurve(20.0)), + fixed = 1.0, + start_up = (hot = 100.0, warm = 250.0, cold = 500.0), + shut_down = 50.0, + ), + base_power = 100.0, + ) + thermal2 = PSY.ThermalMultiStart(; + name = "thermal2", + available = true, + status = true, + bus = bus2, + active_power = 0.5, + reactive_power = 0.0, + rating = 1.0, + prime_mover_type = PSY.PrimeMovers.ST, + fuel = PSY.ThermalFuels.COAL, + active_power_limits = (min = 0.1, max = 0.8), + reactive_power_limits = nothing, + ramp_limits = (up = 0.006, down = 0.006), + power_trajectory = (startup = 0.15, shutdown = 0.15), + time_limits = (up = 2.0, down = 2.0), + start_time_limits = (hot = 1.5, warm = 3.0, cold = 8.0), + start_types = 3, + operation_cost = PSY.ThermalGenerationCost(; + variable = PSY.CostCurve(PSY.QuadraticCurve(5.0, 15.0, 0.0)), + fixed = 0.5, + start_up = (hot = 80.0, warm = 160.0, cold = 320.0), + shut_down = 40.0, + ), + base_power = 100.0, + ) + PSY.add_component!(sys, thermal1) + PSY.add_component!(sys, thermal2) + thermal3 = PSY.ThermalStandard(; + name = "thermal3", + available = true, + status = true, + bus = bus1, + active_power = 0.4, + reactive_power = 0.0, + rating = 0.9, + active_power_limits = (min = 0.15, max = 0.9), + reactive_power_limits = nothing, + ramp_limits = (up = 0.005, down = 0.005), + operation_cost = PSY.ThermalGenerationCost(; + variable = PSY.CostCurve(PSY.LinearCurve(25.0)), + fixed = 0.8, + start_up = 90.0, + shut_down = 45.0, + ), + base_power = 100.0, + time_limits = (up = 2.0, down = 2.0), + ) + PSY.add_component!(sys, thermal3) + renewable = PSY.RenewableDispatch(; + name = "renewable1", + available = true, + bus = bus2, + active_power = 0.0, + reactive_power = 0.0, + rating = 1.2, + prime_mover_type = PSY.PrimeMovers.WT, + reactive_power_limits = nothing, + power_factor = 1.0, + operation_cost = PSY.RenewableGenerationCost(nothing), + base_power = 100.0, + ) + PSY.add_component!(sys, renewable) + load = PSY.PowerLoad(; + name = "load1", + available = true, + bus = bus2, + active_power = 1.0, + reactive_power = 0.0, + base_power = 100.0, + max_active_power = 1.0, + max_reactive_power = 0.0, + ) + PSY.add_component!(sys, load) + + horizon_count = 24 + initial_time = Dates.DateTime("2024-01-01T00:00:00") + load_data = Dict( + initial_time => + [0.6 + 0.4 * abs(sin(pi * t / 12)) for t in 1:horizon_count], + ) + re_data = Dict( + initial_time => + [0.5 + 0.5 * abs(cos(pi * t / 12)) for t in 1:horizon_count], + ) + PSY.add_time_series!( + sys, + load, + PSY.Deterministic(; + name = "max_active_power", + data = load_data, + resolution = Dates.Hour(1), + ), + ) + PSY.add_time_series!( + sys, + renewable, + PSY.Deterministic(; + name = "max_active_power", + data = re_data, + resolution = Dates.Hour(1), + ), + ) + return sys +end + +function _precompile_uc_template() + template = PowerOperationsProblemTemplate(CopperPlateNetworkModel) + set_device_model!(template, PSY.ThermalMultiStart, ThermalMultiStartUnitCommitment) + set_device_model!(template, PSY.ThermalStandard, ThermalStandardUnitCommitment) + set_device_model!(template, PSY.RenewableDispatch, RenewableFullDispatch) + set_device_model!(template, PSY.PowerLoad, StaticPowerLoad) + return template +end + +function _build_precompile_model(sys, template, output_dir) + model = DecisionModel( + template, + sys; + horizon = Dates.Hour(24), + initialize_model = false, + store_variable_names = true, + ) + status = build!( + model; + output_dir = output_dir, + console_level = Logging.Error, + store_system_in_results = false, + ) + if status != IOM.ModelBuildStatus.BUILT + error( + "precompile workload build failed with status $status; " * + "see $(joinpath(output_dir, "operation_problem.log"))", + ) + end + return model +end + +function _precompile_mock_optimizer() + mock = JuMP.MOI.Utilities.MockOptimizer(JuMP.MOI.Utilities.Model{Float64}()) + JuMP.MOI.Utilities.set_mock_optimize!( + mock, + m -> JuMP.MOI.Utilities.mock_optimize!( + m, + JuMP.MOI.OPTIMAL, + ( + JuMP.MOI.FEASIBLE_POINT, + zeros(JuMP.MOI.get(m, JuMP.MOI.NumberOfVariables())), + ), + ), + ) + return mock +end + +function _run_precompile_workload(sys, output_dir) + return Logging.with_logger(Logging.NullLogger()) do + uc_dir = mkpath(joinpath(output_dir, "uc_copperplate")) + model = _build_precompile_model(sys, _precompile_uc_template(), uc_dir) + run_status = solve!( + model; + optimizer = _precompile_mock_optimizer, + console_level = Logging.Error, + export_optimization_problem = false, + ) + if run_status != IOM.RunStatus.SUCCESSFULLY_FINALIZED + error( + "precompile workload solve failed with status $run_status; " * + "see $(joinpath(uc_dir, "operation_problem.log"))", + ) + end + + template_ed_ptdf = PowerOperationsProblemTemplate(PTDFNetworkModel) + set_device_model!(template_ed_ptdf, PSY.ThermalMultiStart, ThermalStandardDispatch) + set_device_model!(template_ed_ptdf, PSY.ThermalStandard, ThermalStandardDispatch) + set_device_model!(template_ed_ptdf, PSY.RenewableDispatch, RenewableFullDispatch) + set_device_model!(template_ed_ptdf, PSY.PowerLoad, StaticPowerLoad) + set_device_model!(template_ed_ptdf, PSY.Line, StaticBranch) + _build_precompile_model( + sys, + template_ed_ptdf, + mkpath(joinpath(output_dir, "ed_ptdf")), + ) + + template_uc_dcp = PowerOperationsProblemTemplate(DCPNetworkModel) + set_device_model!( + template_uc_dcp, + PSY.ThermalMultiStart, + ThermalMultiStartUnitCommitment, + ) + set_device_model!( + template_uc_dcp, + PSY.ThermalStandard, + ThermalStandardUnitCommitment, + ) + set_device_model!(template_uc_dcp, PSY.RenewableDispatch, RenewableFullDispatch) + set_device_model!(template_uc_dcp, PSY.PowerLoad, StaticPowerLoad) + set_device_model!(template_uc_dcp, PSY.Line, StaticBranch) + _build_precompile_model( + sys, + template_uc_dcp, + mkpath(joinpath(output_dir, "uc_dcp")), + ) + + template_ed_areaptdf = PowerOperationsProblemTemplate(AreaPTDFNetworkModel) + set_device_model!( + template_ed_areaptdf, + PSY.ThermalMultiStart, + ThermalStandardDispatch, + ) + set_device_model!( + template_ed_areaptdf, + PSY.ThermalStandard, + ThermalStandardDispatch, + ) + set_device_model!( + template_ed_areaptdf, + PSY.RenewableDispatch, + RenewableFullDispatch, + ) + set_device_model!(template_ed_areaptdf, PSY.PowerLoad, StaticPowerLoad) + set_device_model!(template_ed_areaptdf, PSY.Line, StaticBranchBounds) + _build_precompile_model( + sys, + template_ed_areaptdf, + mkpath(joinpath(output_dir, "ed_areaptdf")), + ) + return nothing + end +end + +PrecompileTools.@setup_workload begin + PrecompileTools.@compile_workload begin + _precompile_output_dir = mktempdir() + try + _precompile_sys = _build_precompile_system() + _run_precompile_workload(_precompile_sys, _precompile_output_dir) + finally + rm(_precompile_output_dir; force = true, recursive = true) + end + end +end diff --git a/test/test_precompile_workload.jl b/test/test_precompile_workload.jl new file mode 100644 index 00000000..a532ac02 --- /dev/null +++ b/test/test_precompile_workload.jl @@ -0,0 +1,49 @@ +@testset "Precompile workload functions run" begin + sys = POM._build_precompile_system() + @test !isempty(PSY.get_components(PSY.ThermalMultiStart, sys)) + @test !isempty(PSY.get_components(PSY.ThermalStandard, sys)) + # AreaPTDFNetworkModel template requires every bus assigned to an Area. + @test length(PSY.get_components(PSY.Area, sys)) == 2 + for bus in PSY.get_components(PSY.ACBus, sys) + @test !isnothing(PSY.get_area(bus)) + end + out = mktempdir(; cleanup = true) + @test POM._run_precompile_workload(sys, out) === nothing +end + +# Ramp and duration constraints are silently dropped for non-binding device +# data (IOM _get_ramp_constraint_devices, POM _get_data_for_tdc). If a future +# edit to the micro system un-binds them, the workload would silently stop +# compiling those paths — this testset makes that a loud failure. +@testset "Precompile workload UC model constructs ramp and duration constraints" begin + sys = POM._build_precompile_system() + model = POM._build_precompile_model( + sys, + POM._precompile_uc_template(), + mktempdir(; cleanup = true), + ) + constraints = IOM.get_constraints(model) + for key in ( + IOM.ConstraintKey(RampConstraint, PSY.ThermalMultiStart, "up"), + IOM.ConstraintKey(RampConstraint, PSY.ThermalMultiStart, "dn"), + IOM.ConstraintKey(DurationConstraint, PSY.ThermalMultiStart, "up"), + IOM.ConstraintKey(DurationConstraint, PSY.ThermalMultiStart, "dn"), + IOM.ConstraintKey(RampConstraint, PSY.ThermalStandard, "up"), + IOM.ConstraintKey(RampConstraint, PSY.ThermalStandard, "dn"), + IOM.ConstraintKey(DurationConstraint, PSY.ThermalStandard, "up"), + IOM.ConstraintKey(DurationConstraint, PSY.ThermalStandard, "dn"), + IOM.ConstraintKey(StartTypeConstraint, PSY.ThermalMultiStart), + IOM.ConstraintKey( + StartupTimeLimitTemperatureConstraint, + PSY.ThermalMultiStart, + "hot", + ), + IOM.ConstraintKey( + StartupTimeLimitTemperatureConstraint, + PSY.ThermalMultiStart, + "warm", + ), + ) + @test haskey(constraints, key) + end +end