diff --git a/LICENSE b/LICENSE index 8b3f7b22d0..f0323d2e2c 100644 --- a/LICENSE +++ b/LICENSE @@ -1,6 +1,6 @@ BSD 3-Clause License -Copyright (c) 2018, 2023 Alliance for Sustainable Energy, LLC +Copyright (c) 2018, 2023 Alliance for Sustainable Energy, LLC and The Regents of the University of California All rights reserved. Redistribution and use in source and binary forms, with or without diff --git a/Project.toml b/Project.toml index 6077b2ca0e..54cbbc0e2a 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "PowerSimulations" uuid = "e690365d-45e2-57bb-ac84-44ba829e73c4" authors = ["Jose Daniel Lara", "Clayton Barrows", "Daniel Thom", "Dheepak Krishnamurthy", "Sourabh Dalvi"] -version = "0.25.2" +version = "0.27.0" [deps] CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b" @@ -46,14 +46,14 @@ JuMP = "1" LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" -PowerModels = "~0.19" -PowerNetworkMatrices = "^0.9" +PowerModels = "^0.20" +PowerNetworkMatrices = "^0.10" PowerFlows = "0.6" PowerSystems = "^3" PrettyTables = "2" ProgressMeter = "^1.5" SHA = "0.7" Serialization = "1" -TimeSeries = "~0.23" +TimeSeries = "~0.23, ~0.24" TimerOutputs = "~0.5" julia = "^1.6" diff --git a/docs/src/tutorials/pcm_simulation.md b/docs/src/tutorials/pcm_simulation.md index 048b49e8dd..bcf7f294fd 100644 --- a/docs/src/tutorials/pcm_simulation.md +++ b/docs/src/tutorials/pcm_simulation.md @@ -78,7 +78,7 @@ example, PSI also provides pre-specified templates for some standard problems: ```@example pcm template_ed = template_economic_dispatch( - network = NetworkModel(StandardPTDFModel, use_slacks = true), + network = NetworkModel(PTDFPowerModel, use_slacks = true), ) ``` diff --git a/src/PowerSimulations.jl b/src/PowerSimulations.jl index f4fb316863..36787698c3 100644 --- a/src/PowerSimulations.jl +++ b/src/PowerSimulations.jl @@ -19,7 +19,6 @@ export SimulationPartitionResults # Network Relevant Exports export NetworkModel -export StandardPTDFModel export PTDFPowerModel export CopperPlatePowerModel export AreaBalancePowerModel diff --git a/src/core/definitions.jl b/src/core/definitions.jl index b51e41e7db..142967f063 100644 --- a/src/core/definitions.jl +++ b/src/core/definitions.jl @@ -53,9 +53,6 @@ const KiB = 1024 const MiB = KiB * KiB const GiB = MiB * KiB -const UNSUPPORTED_POWERMODELS = - [PM.SOCBFPowerModel, PM.SOCBFConicPowerModel, PM.IVRPowerModel] - const PSI_NAME_DELIMITER = "__" const M_VALUE = 1e6 diff --git a/src/core/formulations.jl b/src/core/formulations.jl index 87a6a168c2..b9ea8a7748 100644 --- a/src/core/formulations.jl +++ b/src/core/formulations.jl @@ -152,7 +152,6 @@ struct LossLessLine <: AbstractBranchFormulation end # These formulations are taken directly from PowerModels abstract type AbstractPTDFModel <: PM.AbstractDCPModel end -struct StandardPTDFModel <: AbstractPTDFModel end struct PTDFPowerModel <: AbstractPTDFModel end struct CopperPlatePowerModel <: PM.AbstractActivePowerModel end diff --git a/src/core/network_model.jl b/src/core/network_model.jl index bd165dc63f..6dddb1281e 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -20,11 +20,11 @@ Establishes the model for a particular device specified by type. - `use_slacks::Bool`: Adds slacks to the network modelings - `PTDF::PTDF`: PTDF Array calculated using PowerNetworkMatrices - `duals::Vector{DataType}`: Constraint types to calculate the duals - + - `reduce_radial_branches::Bool`: Skips modeling radial branches in the system to reduce problem size # Example ptdf_array = PTDF(system) -thermal_gens = NetworkModel(StandardPTDFModel, ptdf = ptdf_array), +thermal_gens = NetworkModel(PTDFPowerModel, ptdf = ptdf_array), """ mutable struct NetworkModel{T <: PM.AbstractPowerModel} use_slacks::Bool @@ -32,7 +32,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} subnetworks::Dict{Int, Set{Int}} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} - radial_branches::PNM.RadialBranches + radial_network_reduction::PNM.RadialNetworkReduction reduce_radial_branches::Bool powerflow_evaluation::Union{Nothing, PFS.PowerFlowEvaluationModel} @@ -52,7 +52,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} subnetworks, Dict{PSY.ACBus, Int}(), duals, - PNM.RadialBranches(), + PNM.RadialNetworkReduction(), reduce_radial_branches, powerflow_evaluation, ) @@ -62,7 +62,7 @@ end get_use_slacks(m::NetworkModel) = m.use_slacks get_PTDF_matrix(m::NetworkModel) = m.PTDF_matrix get_reduce_radial_branches(m::NetworkModel) = m.reduce_radial_branches -get_radial_branches(m::NetworkModel) = m.radial_branches +get_radial_network_reduction(m::NetworkModel) = m.radial_network_reduction get_duals(m::NetworkModel) = m.duals get_network_formulation(::NetworkModel{T}) where {T} = T get_reference_buses(m::NetworkModel{T}) where {T <: PM.AbstractPowerModel} = @@ -79,6 +79,15 @@ function add_dual!(model::NetworkModel, dual) return end +function check_radial_branch_reduction_compatibility( + ::Type{T}, +) where {T <: PM.AbstractPowerModel} + if T ∈ INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS + error("Network Model $T is not compatible with radial branch reduction") + end + return +end + function instantiate_network_model( model::NetworkModel{T}, sys::PSY.System, @@ -87,7 +96,8 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches - model.radial_branches = PNM.RadialBranches(sys) + check_radial_branch_reduction_compatibility(T) + model.radial_network_reduction = PNM.RadialNetworkReduction(sys) end return end @@ -107,44 +117,52 @@ function instantiate_network_model( return end -function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys::PSY.System) +function instantiate_network_model(model::NetworkModel{PTDFPowerModel}, sys::PSY.System) if get_PTDF_matrix(model) === nothing @info "PTDF Matrix not provided. Calculating using PowerNetworkMatrices.PTDF" model.PTDF_matrix = PNM.PTDF(sys; reduce_radial_branches = model.reduce_radial_branches) end + if model.reduce_radial_branches + @assert !isempty(model.PTDF_matrix.radial_network_reduction) + model.radial_network_reduction = model.PTDF_matrix.radial_network_reduction + end get_PTDF_matrix(model).subnetworks model.subnetworks = deepcopy(get_PTDF_matrix(model).subnetworks) if length(model.subnetworks) > 1 @debug "System Contains Multiple Subnetworks. Assigning buses to subnetworks." _assign_subnetworks_to_buses(model, sys) end - if model.reduce_radial_branches - @assert !isempty(model.PTDF_matrix.radial_branches) - model.radial_branches = model.PTDF_matrix.radial_branches - end return end function _assign_subnetworks_to_buses( model::NetworkModel{T}, sys::PSY.System, -) where {T <: Union{CopperPlatePowerModel, StandardPTDFModel}} +) where {T <: Union{CopperPlatePowerModel, PTDFPowerModel}} subnetworks = model.subnetworks temp_bus_map = Dict{Int, Int}() - for bus in get_available_components(PSY.ACBus, sys) + radial_network_reduction = PSI.get_radial_network_reduction(model) + for bus in PSI.get_available_components(PSY.ACBus, sys) bus_no = PSY.get_number(bus) + mapped_bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus) if haskey(temp_bus_map, bus_no) model.bus_area_map[bus] = temp_bus_map[bus_no] + continue else + bus_mapped = false for (subnet, bus_set) in subnetworks - if bus_no ∈ bus_set + if mapped_bus_no ∈ bus_set temp_bus_map[bus_no] = subnet model.bus_area_map[bus] = subnet + bus_mapped = true break end end end + if !bus_mapped + error("Bus $(PSY.summary(bus)) not mapped to any reference bus") + end end return end diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 8cb11f374b..a307ed4ab1 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -496,7 +496,7 @@ function _make_system_expressions!( dc_bus_numbers::Vector{Int}, ::Type{T}, bus_reduction_map::Dict{Int64, Set{Int64}}, -) where {T <: Union{PTDFPowerModel, StandardPTDFModel}} +) where {T <: PTDFPowerModel} time_steps = get_time_steps(container) if isempty(bus_reduction_map) ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) @@ -510,6 +510,8 @@ function _make_system_expressions!( ExpressionKey(ActivePowerBalance, PSY.DCBus) => _make_container_array(dc_bus_numbers, time_steps), ExpressionKey(ActivePowerBalance, PSY.ACBus) => + # Bus numbers are sorted to guarantee consistency in the order between the + # containers _make_container_array(sort!(ac_bus_numbers), time_steps), ) return @@ -539,7 +541,7 @@ function build_impl!( transmission, transmission_model.subnetworks, sys, - transmission_model.radial_branches.bus_reduction_map) + transmission_model.radial_network_reduction.bus_reduction_map) # Order is required for device_model in values(template.devices) diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 616500cafb..365d691740 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -211,7 +211,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, device_model::DeviceModel{T, StaticBranch}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(device_model, "filter_function")) @@ -230,7 +230,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, device_model::DeviceModel{T, StaticBranch}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(device_model, "filter_function")) @@ -245,7 +245,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, device_model::DeviceModel{T, StaticBranchBounds}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(device_model, "filter_function")) @@ -264,7 +264,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, device_model::DeviceModel{T, StaticBranchBounds}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(device_model, "filter_function")) @@ -284,7 +284,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, device_model::DeviceModel{T, StaticBranchUnbounded}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(device_model, "filter_function")) @@ -303,7 +303,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, model::DeviceModel{T, StaticBranchUnbounded}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -439,7 +439,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, model::DeviceModel{T, HVDCTwoTerminalUnbounded}, - network_model::NetworkModel{<:Union{StandardPTDFModel, PTDFPowerModel}}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: TwoTerminalHVDCTypes} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -461,7 +461,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, model::DeviceModel{<:TwoTerminalHVDCTypes, HVDCTwoTerminalUnbounded}, - network_model::NetworkModel{<:Union{StandardPTDFModel, PTDFPowerModel}}, + network_model::NetworkModel{PTDFPowerModel}, ) add_constraint_dual!(container, sys, model) return @@ -497,7 +497,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, model::DeviceModel{T, HVDCTwoTerminalLossless}, - network_model::NetworkModel{<:Union{StandardPTDFModel, PTDFPowerModel}}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: TwoTerminalHVDCTypes} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -522,7 +522,7 @@ function construct_device!( network_model::NetworkModel{U}, ) where { T <: TwoTerminalHVDCTypes, - U <: Union{StandardPTDFModel, PTDFPowerModel}, + U <: PTDFPowerModel, } devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -536,7 +536,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, model::DeviceModel{T, HVDCTwoTerminalDispatch}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: TwoTerminalHVDCTypes} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -586,7 +586,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, model::DeviceModel{T, HVDCTwoTerminalDispatch}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: TwoTerminalHVDCTypes} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) @@ -705,7 +705,7 @@ function construct_device!( sys::PSY.System, ::ArgumentConstructStage, model::DeviceModel{PSY.PhaseShiftingTransformer, PhaseAngleControl}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) devices = get_available_components( PSY.PhaseShiftingTransformer, @@ -749,7 +749,7 @@ function construct_device!( sys::PSY.System, ::ModelConstructStage, model::DeviceModel{PSY.PhaseShiftingTransformer, PhaseAngleControl}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) devices = get_available_components( PSY.PhaseShiftingTransformer, diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index de74ab915c..64d5e47f54 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -46,7 +46,7 @@ end function add_variables!( container::OptimizationContainer, ::Type{FlowActivePowerVariable}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, devices::IS.FlattenIteratorWrapper{T}, formulation::AbstractBranchFormulation, ) where {T <: PSY.ACBranch} @@ -108,8 +108,8 @@ function branch_rate_bounds!( ) where {B <: PSY.ACBranch} var = get_variable(container, FlowActivePowerVariable(), B) - radial_branches = get_radial_branches(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices name = PSY.get_name(d) @@ -136,8 +136,8 @@ function branch_rate_bounds!( ] time_steps = get_time_steps(container) - radial_branches = get_radial_branches(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices name = PSY.get_name(d) @@ -191,11 +191,11 @@ function add_constraints!( V <: PM.AbstractActivePowerModel, } time_steps = get_time_steps(container) - radial_branches = get_radial_branches(network_model) - if isempty(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + if isempty(radial_network_reduction) device_names = [PSY.get_name(d) for d in devices] else - device_names = PNM.get_meshed_branches(radial_branches) + device_names = PNM.get_meshed_branches(radial_network_reduction) end con_lb = @@ -221,7 +221,7 @@ function add_constraints!( for device in devices ci_name = PSY.get_name(device) - if ci_name ∈ PNM.get_radial_branches(radial_branches) + if ci_name ∈ PNM.get_radial_branches(radial_network_reduction) continue end limits = get_min_max_limits(device, RateLimitConstraint, U) # depends on constraint type and formulation type @@ -287,8 +287,8 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_branches = get_radial_branches(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data if r[1] ∈ radial_branches_names @@ -328,8 +328,8 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_branches = get_radial_branches(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data if r[1] ∈ radial_branches_names @@ -345,14 +345,14 @@ function add_constraints!( end """ -Add network flow constraints for ACBranch and NetworkModel with StandardPTDFModel +Add network flow constraints for ACBranch and NetworkModel with PTDFPowerModel """ function add_constraints!( container::OptimizationContainer, ::Type{NetworkFlowConstraint}, devices::IS.FlattenIteratorWrapper{B}, model::DeviceModel{B, <:AbstractBranchFormulation}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {B <: PSY.ACBranch} ptdf = get_PTDF_matrix(network_model) # This is a workaround to not call the same list comprehension to find @@ -388,14 +388,14 @@ function add_constraints!( end """ -Add network flow constraints for PhaseShiftingTransformer and NetworkModel with StandardPTDFModel +Add network flow constraints for PhaseShiftingTransformer and NetworkModel with PTDFPowerModel """ function add_constraints!( container::OptimizationContainer, ::Type{NetworkFlowConstraint}, devices::IS.FlattenIteratorWrapper{T}, model::DeviceModel{T, PhaseAngleControl}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: PSY.PhaseShiftingTransformer} ptdf = get_PTDF_matrix(network_model) branches = PSY.get_name.(devices) diff --git a/src/devices_models/devices/HVDCsystems.jl b/src/devices_models/devices/HVDCsystems.jl index ee20e478f6..ad36f26954 100644 --- a/src/devices_models/devices/HVDCsystems.jl +++ b/src/devices_models/devices/HVDCsystems.jl @@ -215,7 +215,7 @@ function add_to_expression!( U <: ActivePowerVariable, V <: PSY.InterconnectingConverter, W <: AbstractConverterFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel } variable = get_variable(container, U(), V) expression_dc = get_expression(container, T(), PSY.DCBus) diff --git a/src/devices_models/devices/common/add_constraint_dual.jl b/src/devices_models/devices/common/add_constraint_dual.jl index 1d5e9111bd..0b49570ba2 100644 --- a/src/devices_models/devices/common/add_constraint_dual.jl +++ b/src/devices_models/devices/common/add_constraint_dual.jl @@ -31,7 +31,7 @@ function add_constraint_dual!( container::OptimizationContainer, sys::PSY.System, model::NetworkModel{T}, -) where {T <: Union{CopperPlatePowerModel, StandardPTDFModel}} +) where {T <: Union{CopperPlatePowerModel, PTDFPowerModel}} if !isempty(get_duals(model)) for constraint_type in get_duals(model) assign_dual_variable!(container, constraint_type, sys, model) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 17bdefdbb0..29bd067dd9 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -89,9 +89,9 @@ function add_to_expression!( } param_container = get_parameter(container, U(), V) multiplier = get_multiplier_array(param_container) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices, t in get_time_steps(container) - bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) name = PSY.get_name(d) _add_to_jump_expression!( get_expression(container, T(), PSY.ACBus)[bus_no, t], @@ -117,9 +117,9 @@ function add_to_expression!( X <: PM.AbstractPowerModel, } parameter = get_parameter_array(container, U(), V) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices, t in get_time_steps(container) - bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) name = PSY.get_name(d) mult = get_expression_multiplier(U(), T(), d, W()) _add_to_jump_expression!( @@ -150,10 +150,10 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) - bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) _add_to_jump_expression!( expression[bus_no, t], variable[name, t], @@ -178,7 +178,7 @@ function add_to_expression!( U <: HVDCLosses, V <: TwoTerminalHVDCTypes, W <: HVDCTwoTerminalDispatch, - X <: Union{StandardPTDFModel, CopperPlatePowerModel}, + X <: Union{PTDFPowerModel, CopperPlatePowerModel}, } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.System) @@ -212,15 +212,14 @@ function add_to_expression!( U <: FlowActivePowerToFromVariable, V <: TwoTerminalHVDCTypes, W <: AbstractDeviceFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } var = get_variable(container, U(), V) nodal_expr = get_expression(container, T(), PSY.ACBus) sys_expr = get_expression(container, T(), PSY.System) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_to_) + bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).to) ref_bus_from = get_reference_bus(network_model, PSY.get_arc(d).from) ref_bus_to = get_reference_bus(network_model, PSY.get_arc(d).to) for t in get_time_steps(container) @@ -249,15 +248,15 @@ function add_to_expression!( U <: FlowActivePowerFromToVariable, V <: TwoTerminalHVDCTypes, W <: AbstractTwoTerminalDCLineFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } var = get_variable(container, U(), V) nodal_expr = get_expression(container, T(), PSY.ACBus) sys_expr = get_expression(container, T(), PSY.System) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) - bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_from_) + bus_no_from = + PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).from) ref_bus_to = get_reference_bus(network_model, PSY.get_arc(d).to) ref_bus_from = get_reference_bus(network_model, PSY.get_arc(d).from) for t in get_time_steps(container) @@ -358,11 +357,11 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_arc(d).from) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( expression[bus_no, t], @@ -393,11 +392,11 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_arc(d).to) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( expression[bus_no, t], @@ -425,11 +424,11 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) _add_to_jump_expression!( expression[bus_no, t], variable[name, t], @@ -580,18 +579,18 @@ function add_to_expression!( U <: TimeSeriesParameter, V <: PSY.StaticInjection, W <: AbstractDeviceFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } param_container = get_parameter(container, U(), V) multiplier = get_multiplier_array(param_container) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) bus_no_ = PSY.get_number(device_bus) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) ref_bus = get_reference_bus(network_model, device_bus) param = get_parameter_column_refs(param_container, name) for t in get_time_steps(container) @@ -614,16 +613,16 @@ function add_to_expression!( U <: OnStatusParameter, V <: PSY.ThermalGen, W <: AbstractDeviceFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } parameter = get_parameter_array(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) mult = get_expression_multiplier(U(), T(), d, W()) device_bus = PSY.get_bus(d) ref_bus = get_reference_bus(network_model, device_bus) @@ -648,17 +647,16 @@ function add_to_expression!( U <: VariableType, V <: PSY.StaticInjection, W <: AbstractDeviceFormulation, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } variable = get_variable(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no_ = PSY.get_number(device_bus) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, device_bus) ref_bus = get_reference_bus(network_model, device_bus) for t in get_time_steps(container) _add_to_jump_expression!( @@ -688,18 +686,16 @@ function add_to_expression!( U <: OnVariable, V <: PSY.ThermalGen, W <: Union{AbstractCompactUnitCommitment, ThermalCompactDispatch}, - X <: Union{PTDFPowerModel, StandardPTDFModel}, + X <: PTDFPowerModel, } variable = get_variable(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices name = PSY.get_name(d) - device_bus = PSY.get_bus(d) - bus_no_ = PSY.get_number(device_bus) - bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) - ref_bus = get_reference_bus(network_model, device_bus) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) + ref_bus = get_reference_bus(network_model, PSY.get_bus(d)) for t in get_time_steps(container) _add_to_jump_expression!( sys_expr[ref_bus, t], @@ -735,12 +731,11 @@ function add_to_expression!( } var = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) - bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_to_) + bus_no_from = + PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).from) + bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).to) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( @@ -767,7 +762,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where { T <: ActivePowerBalance, U <: FlowActivePowerVariable, @@ -777,12 +772,11 @@ function add_to_expression!( var = get_variable(container, U(), V) nodal_expr = get_expression(container, T(), PSY.ACBus) sys_expr = get_expression(container, T(), PSY.System) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) - bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_to_) + bus_no_from = + PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).from) + bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).to) ref_bus_from = get_reference_bus(network_model, PSY.get_arc(d).from) ref_bus_to = get_reference_bus(network_model, PSY.get_arc(d).to) for t in get_time_steps(container) @@ -850,16 +844,15 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{PSY.PhaseShiftingTransformer}, ::DeviceModel{PSY.PhaseShiftingTransformer, V}, - network_model::NetworkModel{StandardPTDFModel}, + network_model::NetworkModel{PTDFPowerModel}, ) where {T <: ActivePowerBalance, U <: PhaseShifterAngle, V <: PhaseAngleControl} var = get_variable(container, U(), PSY.PhaseShiftingTransformer) expression = get_expression(container, T(), PSY.ACBus) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) - bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_to_) + bus_no_from = + PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).from) + bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).to) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( @@ -1083,7 +1076,7 @@ function add_to_expression!( ) where { T <: ActivePowerBalance, U <: Union{SystemBalanceSlackUp, SystemBalanceSlackDown}, - W <: Union{CopperPlatePowerModel, StandardPTDFModel}, + W <: Union{CopperPlatePowerModel, PTDFPowerModel}, } variable = get_variable(container, U(), PSY.System) expression = get_expression(container, T(), PSY.System) @@ -1111,8 +1104,7 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus) expression = get_expression(container, T(), PSY.ACBus) - #TODO: Update Slacks to consider reduced buses - #@assert_op length(axes(variable, 1)) == length(axes(expression, 1)) + @assert_op length(axes(variable, 1)) == length(axes(expression, 1)) # We uses axis here to avoid double addition of the slacks to the aggregated buses for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index b58d01f965..95e3733aac 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -341,11 +341,21 @@ end Startup and shutdown active power limits for Compact Unit Commitment """ function get_startup_shutdown_limits( - device, + device::PSY.ThermalMultiStart, ::Type{ActivePowerVariableLimitsConstraint}, - ::Type{<:ThermalMultiStartUnitCommitment}, + ::Type{ThermalMultiStartUnitCommitment}, ) - return PSY.get_power_trajectory(device) + startup_shutdown = PSY.get_power_trajectory(device) + if isnothing(startup_shutdown) + @warn( + "Generator $(summary(device)) has a Nothing startup_shutdown property. Using active power limits." + ) + return ( + startup = PSY.get_active_power_limits(device).max, + shutdown = PSY.get_active_power_limits(device).max, + ) + end + return startup_shutdown end """ @@ -355,7 +365,7 @@ function get_min_max_limits( device, ::Type{ActivePowerVariableLimitsConstraint}, ::Type{<:AbstractCompactUnitCommitment}, -) # -> Union{Nothing, NamedTuple{(:startup, :shutdown), Tuple{Float64, Float64}}} +) # -> Union{Nothing, NamedTuple{(:min, :max), Tuple{Float64, Float64}}} return ( min = 0, max = PSY.get_active_power_limits(device).max - @@ -562,6 +572,7 @@ function add_constraints!( name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type startup_shutdown_limits = get_startup_shutdown_limits(device, T, W) + @assert !isnothing(startup_shutdown_limits) "$(name)" if JuMP.has_lower_bound(varp[name, t]) JuMP.set_lower_bound(varp[name, t], 0.0) end diff --git a/src/initial_conditions/initialization.jl b/src/initial_conditions/initialization.jl index a2158881f3..01ec41fbec 100644 --- a/src/initial_conditions/initialization.jl +++ b/src/initial_conditions/initialization.jl @@ -9,7 +9,8 @@ function get_initial_conditions_template(model::OperationModel) get_network_model(model.template), ), ) - network_model.radial_branches = get_radial_branches(get_network_model(model.template)) + network_model.radial_network_reduction = + get_radial_network_reduction(get_network_model(model.template)) network_model.subnetworks = get_subnetworks(get_network_model(model.template)) # Initialization does not support PowerFlow evaluation network_model.powerflow_evaluation = nothing diff --git a/src/network_models/copperplate_model.jl b/src/network_models/copperplate_model.jl index 3d4c12b01b..1cc0886148 100644 --- a/src/network_models/copperplate_model.jl +++ b/src/network_models/copperplate_model.jl @@ -6,7 +6,7 @@ function add_constraints!( ) where { T <: CopperPlateBalanceConstraint, U <: PSY.System, - V <: Union{CopperPlatePowerModel, StandardPTDFModel, PTDFPowerModel}, + V <: Union{CopperPlatePowerModel, PTDFPowerModel}, } time_steps = get_time_steps(container) expressions = get_expression(container, ActivePowerBalance(), U) diff --git a/src/network_models/network_constructor.jl b/src/network_models/network_constructor.jl index 7f53454439..890c6d98d9 100644 --- a/src/network_models/network_constructor.jl +++ b/src/network_models/network_constructor.jl @@ -53,7 +53,7 @@ end function construct_network!( container::OptimizationContainer, sys::PSY.System, - model::NetworkModel{StandardPTDFModel}, + model::NetworkModel{PTDFPowerModel}, ::ProblemTemplate, ) if get_use_slacks(model) @@ -75,34 +75,11 @@ function construct_network!( return end -function construct_network!( - container::OptimizationContainer, - sys::PSY.System, - model::NetworkModel{T}, - template::ProblemTemplate, -) where {T <: PTDFPowerModel} - construct_network!( - container, - sys, - model, - template; - instantiate_model = instantiate_nip_ptdf_expr_model, - ) - - add_pm_expr_refs!(container, T, sys) - - add_constraints!(container, CopperPlateBalanceConstraint, sys, model) - add_constraints!(container, NodalBalanceActiveConstraint, sys, model) - add_constraint_dual!(container, sys, model) - return -end - function construct_network!( container::OptimizationContainer, sys::PSY.System, model::NetworkModel{T}, template::ProblemTemplate; - instantiate_model = instantiate_nip_expr_model, ) where {T <: PM.AbstractActivePowerModel} if T in UNSUPPORTED_POWERMODELS throw( @@ -126,9 +103,9 @@ function construct_network!( objective_function!(container, PSY.ACBus, model) end - @debug "Building the $T network with $instantiate_model method" _group = + @debug "Building the $T network with instantiate_nip_expr_model method" _group = LOG_GROUP_NETWORK_CONSTRUCTION - powermodels_network!(container, T, sys, template, instantiate_model) + powermodels_network!(container, T, sys, template, instantiate_nip_expr_model) #Constraints in case the model has DC Buses add_constraints!(container, NodalBalanceActiveConstraint, sys, model) add_pm_variable_refs!(container, T, sys) @@ -143,7 +120,6 @@ function construct_network!( sys::PSY.System, model::NetworkModel{T}, template::ProblemTemplate; - instantiate_model = instantiate_nip_expr_model, ) where {T <: PM.AbstractPowerModel} if T in UNSUPPORTED_POWERMODELS throw( @@ -174,11 +150,11 @@ function construct_network!( objective_function!(container, PSY.ACBus, model) end - @debug "Building the $T network with $instantiate_model method" _group = + @debug "Building the $T network with instantiate_nip_expr_model method" _group = LOG_GROUP_NETWORK_CONSTRUCTION #Constraints in case the model has DC Buses add_constraints!(container, NodalBalanceActiveConstraint, sys, model) - powermodels_network!(container, T, sys, template, instantiate_model) + powermodels_network!(container, T, sys, template, instantiate_nip_expr_model) add_pm_variable_refs!(container, T, sys) add_pm_constraint_refs!(container, T, sys) @@ -190,8 +166,7 @@ function construct_network!( container::OptimizationContainer, sys::PSY.System, model::NetworkModel{T}, - template::ProblemTemplate; - instantiate_model = instantiate_bfp_expr_model, + template::ProblemTemplate, ) where {T <: PM.AbstractBFModel} if T in UNSUPPORTED_POWERMODELS throw( @@ -235,23 +210,24 @@ function construct_network!( objective_function!(container, PSY.ACBus, model) end - @debug "Building the $T network with $instantiate_model method" _group = + @debug "Building the $T network with instantiate_bfp_expr_model method" _group = LOG_GROUP_NETWORK_CONSTRUCTION #Constraints in case the model has DC Buses add_constraints!(container, NodalBalanceActiveConstraint, sys, model) - powermodels_network!(container, T, sys, template, instantiate_model) + powermodels_network!(container, T, sys, template, instantiate_bfp_expr_model) add_pm_variable_refs!(container, T, sys) add_pm_constraint_refs!(container, T, sys) add_constraint_dual!(container, sys, model) return end +#= +# AbstractIVRModel models not currently supported function construct_network!( container::OptimizationContainer, sys::PSY.System, model::NetworkModel{T}, template::ProblemTemplate; - instantiate_model = instantiate_vip_expr_model, ) where {T <: PM.AbstractIVRModel} if T in UNSUPPORTED_POWERMODELS throw( @@ -299,13 +275,14 @@ function construct_network!( objective_function!(container, PSY.ACBus, model) end - @debug "Building the $T network with $instantiate_model method" _group = + @debug "Building the $T network with instantiate_vip_expr_model method" _group = LOG_GROUP_NETWORK_CONSTRUCTION #Constraints in case the model has DC Buses add_constraints!(container, NodalBalanceActiveConstraint, sys, model) - powermodels_network!(container, T, sys, template, instantiate_model) + powermodels_network!(container, T, sys, template, instantiate_vip_expr_model) add_pm_variable_refs!(container, T, sys) add_pm_constraint_refs!(container, T, sys) add_constraint_dual!(container, sys, model) return end +=# diff --git a/src/network_models/network_slack_variables.jl b/src/network_models/network_slack_variables.jl index 619aeab047..8cc2f16844 100644 --- a/src/network_models/network_slack_variables.jl +++ b/src/network_models/network_slack_variables.jl @@ -10,7 +10,7 @@ function add_variables!( network_model::NetworkModel{U}, ) where { T <: Union{SystemBalanceSlackUp, SystemBalanceSlackDown}, - U <: Union{CopperPlatePowerModel, StandardPTDFModel}, + U <: Union{CopperPlatePowerModel, PTDFPowerModel}, } time_steps = get_time_steps(container) reference_buses = get_reference_buses(network_model) @@ -37,7 +37,13 @@ function add_variables!( U <: PM.AbstractActivePowerModel, } time_steps = get_time_steps(container) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + radial_network_reduction = get_radial_network_reduction(network_model) + if isempty(radial_network_reduction) + bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_numbers = collect(keys(PNM.get_bus_reduction_map(radial_network_reduction))) + end + variable = add_variable_container!(container, T(), PSY.ACBus, bus_numbers, time_steps) for t in time_steps, n in bus_numbers variable[n, t] = JuMP.@variable( @@ -59,7 +65,12 @@ function add_variables!( U <: PM.AbstractPowerModel, } time_steps = get_time_steps(container) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + radial_network_reduction = get_radial_network_reduction(network_model) + if isempty(radial_network_reduction) + bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_numbers = collect(keys(PNM.get_bus_reduction_map(radial_network_reduction))) + end variable_active = add_variable_container!(container, T(), PSY.ACBus, "P", bus_numbers, time_steps) variable_reactive = @@ -84,7 +95,7 @@ function objective_function!( container::OptimizationContainer, ::Type{PSY.System}, network_model::NetworkModel{T}, -) where {T <: Union{CopperPlatePowerModel, StandardPTDFModel}} +) where {T <: Union{CopperPlatePowerModel, PTDFPowerModel}} variable_up = get_variable(container, SystemBalanceSlackUp(), PSY.System) variable_dn = get_variable(container, SystemBalanceSlackDown(), PSY.System) reference_buses = get_reference_buses(network_model) diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index ec012ff097..e2dc9cdac7 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -389,8 +389,8 @@ function get_branches_to_pm( PM_branches = Dict{String, Any}() PMmap_br = Dict{PM_MAP_TUPLE, T}() - radial_branches = get_radial_branches(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for (d, device_model) in branch_template comp_type = get_component_type(device_model) if comp_type <: TwoTerminalHVDCTypes @@ -401,7 +401,7 @@ function get_branches_to_pm( filter_func = get_attribute(device_model, "filter_function") for (i, branch) in enumerate(get_available_components(comp_type, sys, filter_func)) if PSY.get_name(branch) ∈ radial_branches_names - @debug "Skipping branch $(PSY.get_name(branch)) since its radial" + @debug "Skipping branch $(PSY.get_name(branch)) since it is radial" continue end ix = i + start_idx @@ -446,18 +446,6 @@ function get_branches_to_pm( return PM_branches, PMmap_br end -function get_branches_to_pm( - ::PSY.System, - network_model::NetworkModel{PTDFPowerModel}, - ::Type{T}, - branch_template::BranchModelContainer, - start_idx = 0, -) where {T <: TwoTerminalHVDCTypes} - PM_branches = Dict{String, Any}() - PMmap_br = Dict{PM_MAP_TUPLE, T}() - return PM_branches, PMmap_br -end - function get_buses_to_pm(buses::IS.FlattenIteratorWrapper{PSY.ACBus}) PM_buses = Dict{String, Any}() PMmap_buses = Dict{Int, PSY.ACBus}() diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index 942a735985..0e659801ca 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -6,6 +6,13 @@ ################################################################################# # Model Definitions +const UNSUPPORTED_POWERMODELS = + [PM.SOCBFPowerModel, PM.SOCBFConicPowerModel, PM.IVRPowerModel] + +const INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS = [ + PM.SparseSDPWRMPowerModel, +] + function instantiate_nip_expr_model(data::Dict{String, Any}, model_constructor; kwargs...) return PM.instantiate_model(data, model_constructor, instantiate_nip_expr; kwargs...) end @@ -13,7 +20,6 @@ end # replicates PM.build_mn_opf function instantiate_nip_expr(pm::PM.AbstractPowerModel) for n in eachindex(PM.nws(pm)) - @assert !PM.ismulticonductor(pm; nw = n) PM.variable_bus_voltage(pm; nw = n) PM.variable_branch_power(pm; nw = n, bounded = false) PM.variable_dcline_power(pm; nw = n, bounded = false) @@ -43,60 +49,6 @@ function instantiate_nip_expr(pm::PM.AbstractPowerModel) return end -function instantiate_nip_ptdf_expr_model( - data::Dict{String, Any}, - model_constructor; - kwargs..., -) - return PM.instantiate_model( - data, - PM.DCPPowerModel, - instantiate_nip_ptdf_expr; - ref_extensions = [PM.ref_add_connected_components!, PM.ref_add_sm!], - kwargs..., - ) -end - -# replicates PM.build_opf_ptdf -function instantiate_nip_ptdf_expr(pm::PM.AbstractPowerModel) - for n in eachindex(PM.nws(pm)) - @assert !PM.ismulticonductor(pm; nw = n) - - #PM.variable_gen_power(pm) #connect P__* with these - - for i in PM.ids(pm, :bus; nw = n) - if !haskey(PM.var(pm, n), :inj_p) - PM.var(pm, n)[:inj_p] = Dict{Int, Any}() - end - PM.var(pm, n, :inj_p)[i] = PM.ref(pm, :bus, i; nw = n)["inj_p"] # use :nodal_balance_expr - end - - PM.constraint_model_voltage(pm; nw = n) - - # this constraint is implicit in this model - #for i in PM.ids(pm, :ref_buses, nw = n) - # PM.constraint_theta_ref(pm, i, nw = n) - #end - - # done later with PSI copper_plate - #for i in PM.ids(pm, :bus, nw = n) - # constraint_power_balance_ni_expr(pm, i, nw = n) - #end - for (i, branch) in PM.ref(pm, :branch; nw = n) - if haskey(branch, "rate_a") - PM.expression_branch_power_ohms_yt_from_ptdf(pm, i; nw = n) - PM.expression_branch_power_ohms_yt_to_ptdf(pm, i; nw = n) - end - - # done in PSI construct_branch! - #PM.constraint_thermal_limit_from(pm, i, nw = n) - #PM.constraint_thermal_limit_to(pm, i, nw = n) - end - end - - return -end - function instantiate_bfp_expr_model(data::Dict{String, Any}, model_constructor; kwargs...) return PM.instantiate_model(data, model_constructor, instantiate_bfp_expr; kwargs...) end @@ -104,7 +56,6 @@ end # replicates PM.build_mn_opf_bf_strg function instantiate_bfp_expr(pm::PM.AbstractPowerModel) for n in eachindex(PM.nws(pm)) - @assert !PM.ismulticonductor(pm; nw = n) PM.variable_bus_voltage(pm; nw = n) PM.variable_branch_power(pm; nw = n, bounded = false) PM.variable_dcline_power(pm; nw = n, bounded = false) @@ -134,9 +85,12 @@ function instantiate_bfp_expr(pm::PM.AbstractPowerModel) return end +#= +# VI Methdos not supported currently function instantiate_vip_expr_model(data::Dict{String, Any}, model_constructor; kwargs...) throw(error("VI Models not currently supported")) end +=# ################################################################################# # Model Extention Functions @@ -198,6 +152,8 @@ function constraint_power_balance_ni_expr( return end +#= +# VI Methdos not supported currently function constraint_current_balance_ni_expr( pm::PM.AbstractPowerModel, i::Int; @@ -254,6 +210,7 @@ function constraint_current_balance_ni_expr( return end +=# """ active power only models ignore reactive power variables @@ -287,17 +244,17 @@ function powermodels_network!( system_formulation::Type{S}, sys::PSY.System, template::ProblemTemplate, - instantiate_model = instantiate_nip_expr_model, + instantiate_model, ) where {S <: PM.AbstractPowerModel} time_steps = get_time_steps(container) pm_data, PM_map = pass_to_pm(sys, template, time_steps[end]) network_model = get_network_model(template) - radial_branches = get_radial_branches(network_model) - if isempty(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + if isempty(radial_network_reduction) ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) else - bus_reduction_map = radial_branches.bus_reduction_map + bus_reduction_map = PNM.get_bus_reduction_map(radial_network_reduction) ac_bus_numbers = collect(keys(bus_reduction_map)) end @@ -326,18 +283,18 @@ function powermodels_network!( system_formulation::Type{S}, sys::PSY.System, template::ProblemTemplate, - instantiate_model = instantiate_nip_expr_model, + instantiate_model, ) where {S <: PM.AbstractActivePowerModel} time_steps = get_time_steps(container) pm_data, PM_map = pass_to_pm(sys, template, time_steps[end]) buses = get_available_components(PSY.ACBus, sys) network_model = get_network_model(template) - radial_branches = get_radial_branches(network_model) - if isempty(radial_branches) + radial_network_reduction = get_radial_network_reduction(network_model) + if isempty(radial_network_reduction) ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) else - bus_reduction_map = radial_branches.bus_reduction_map + bus_reduction_map = PNM.get_bus_reduction_map(radial_network_reduction) ac_bus_numbers = collect(keys(bus_reduction_map)) end @@ -390,16 +347,6 @@ function PMvarmap(::Type{S}) where {S <: PM.AbstractActivePowerModel} return pm_variable_map end -function PMvarmap(::Type{PTDFPowerModel}) - pm_variable_map = Dict{Type, Dict{Symbol, Union{String, NamedTuple}}}() - - pm_variable_map[PSY.ACBus] = Dict() - pm_variable_map[PSY.ACBranch] = Dict() - pm_variable_map[TwoTerminalHVDCTypes] = Dict() - - return pm_variable_map -end - function PMvarmap(::Type{S}) where {S <: PM.AbstractPowerModel} pm_variable_map = Dict{Type, Dict{Symbol, Union{VariableType, NamedTuple}}}() @@ -454,23 +401,6 @@ function PMexprmap(::Type{S}) where {S <: PM.AbstractPowerModel} return pm_expr_map end -function PMexprmap(::Type{PTDFPowerModel}) - pm_expr_map = Dict{ - Type, - NamedTuple{ - (:pm_expr, :psi_con), - Tuple{Dict{Symbol, Union{VariableType, NamedTuple}}, ConstraintType}, - }, - }() - - pm_expr_map[PSY.ACBranch] = ( - pm_expr = Dict(:p => (from_to = FlowActivePowerVariable(), to_from = nothing)), - psi_con = NetworkFlowConstraint(), - ) - - return pm_expr_map -end - function add_pm_variable_refs!( container::OptimizationContainer, system_formulation::Type{S}, @@ -581,91 +511,3 @@ function add_pm_constraint_refs!( end end end - -function add_pm_expr_refs!( - container::OptimizationContainer, - system_formulation::Type{S}, - ::PSY.System, -) where {S <: PM.AbstractPowerModel} - time_steps = get_time_steps(container) - - ACbranch_dict = container.pm.ext[:PMmap].arcs - ACbranch_types = typeof.(values(ACbranch_dict)) - - pm_variable_types = keys(PM.var(container.pm, 1)) - pm_expr_map = PMexprmap(system_formulation) - - add_pm_expr_refs!( - container, - PSY.ACBranch, - ACbranch_types, - ACbranch_dict, - pm_expr_map, - pm_variable_types, - time_steps, - ) - return -end - -function add_pm_expr_refs!( - container::OptimizationContainer, - d_class::Type, - device_types::Vector, - pm_map::Dict, - pm_expr_map::Dict, - pm_variable_types::Base.KeySet, - time_steps::UnitRange{Int}, -) - for d_type in Set(device_types) - for (pm_expr_var, ps_v) in pm_expr_map[d_class].pm_expr - if pm_expr_var in pm_variable_types - pm_devices = keys(PM.var(container.pm, pm_expr_var; nw = 1)) - mapped_pm_devices = Vector() - mapped_ps_devices = Vector{d_type}() - for d in pm_map - if typeof(d[2]) == d_type && - d[1].from_to ∈ pm_devices && - d[1].to_from ∈ pm_devices - push!(mapped_pm_devices, d[1]) - push!(mapped_ps_devices, d[2]) - end - end - isempty(mapped_ps_devices) && continue - mapped_ps_device_names = PSY.get_name.(mapped_ps_devices) - - # add variable in psi - # add psi_var = pm_expr_var as constraint - for dir in fieldnames(typeof(ps_v)) - var_type = getfield(ps_v, dir) - var_type === nothing && continue - - add_variable!( - container, - var_type, - mapped_ps_devices, - StaticBranchUnbounded(), - ) - psi_variable_container = get_variable(container, var_type, d_type) - - con_type = pm_expr_map[d_class].psi_con - psi_constraint_container = add_constraints_container!( - container, - con_type, - d_type, - mapped_ps_device_names, - time_steps, - ) - for t in time_steps, - (pm_d, name) in zip(mapped_pm_devices, mapped_ps_device_names) - - pm_expr = PM.var(container.pm, t, pm_expr_var, getfield(pm_d, dir)) - psi_constraint_container[name, t] = JuMP.@constraint( - container.JuMPmodel, - psi_variable_container[name, t] == pm_expr - ) - end - end - end - end - end -end diff --git a/src/operation/decision_model.jl b/src/operation/decision_model.jl index 8287346b53..a0cde57e3d 100644 --- a/src/operation/decision_model.jl +++ b/src/operation/decision_model.jl @@ -9,6 +9,15 @@ Generic PowerSimulations Operation Problem Type for unspecified models """ struct GenericOpProblem <: DefaultDecisionProblem end +mutable struct DecisionModel{M <: DecisionProblem} <: OperationModel + name::Symbol + template::ProblemTemplate + sys::PSY.System + internal::Union{Nothing, ModelInternal} + store::DecisionModelStore + ext::Dict{String, Any} +end + """ DecisionModel{M}( template::ProblemTemplate, @@ -51,40 +60,31 @@ template = ProblemTemplate(CopperPlatePowerModel, devices, branches, services) OpModel = DecisionModel(MockOperationProblem, template, system) ``` """ -mutable struct DecisionModel{M <: DecisionProblem} <: OperationModel - name::Symbol - template::ProblemTemplate - sys::PSY.System - internal::Union{Nothing, ModelInternal} - store::DecisionModelStore - ext::Dict{String, Any} - - function DecisionModel{M}( - template::ProblemTemplate, - sys::PSY.System, - settings::Settings, - jump_model::Union{Nothing, JuMP.Model} = nothing; - name = nothing, - ) where {M <: DecisionProblem} - if name === nothing - name = nameof(M) - elseif name isa String - name = Symbol(name) - end - internal = ModelInternal( - OptimizationContainer(sys, settings, jump_model, PSY.Deterministic), - ) - template_ = deepcopy(template) - finalize_template!(template_, sys) - return new{M}( - name, - template_, - sys, - internal, - DecisionModelStore(), - Dict{String, Any}(), - ) +function DecisionModel{M}( + template::ProblemTemplate, + sys::PSY.System, + settings::Settings, + jump_model::Union{Nothing, JuMP.Model} = nothing; + name = nothing, +) where {M <: DecisionProblem} + if name === nothing + name = nameof(M) + elseif name isa String + name = Symbol(name) end + internal = ModelInternal( + OptimizationContainer(sys, settings, jump_model, PSY.Deterministic), + ) + template_ = deepcopy(template) + finalize_template!(template_, sys) + return DecisionModel{M}( + name, + template_, + sys, + internal, + DecisionModelStore(), + Dict{String, Any}(), + ) end function DecisionModel{M}( diff --git a/src/operation/decision_model_store.jl b/src/operation/decision_model_store.jl index b8c2591164..686a5b9586 100644 --- a/src/operation/decision_model_store.jl +++ b/src/operation/decision_model_store.jl @@ -27,7 +27,7 @@ end function initialize_storage!( store::DecisionModelStore, - container::OptimizationContainer, + container::AbstractModelContainer, params::ModelStoreParams, ) num_of_executions = get_num_executions(params) diff --git a/src/operation/model_internal.jl b/src/operation/model_internal.jl index f00b690259..4753e6f648 100644 --- a/src/operation/model_internal.jl +++ b/src/operation/model_internal.jl @@ -12,9 +12,9 @@ mutable struct SimulationInfo sequence_uuid::Base.UUID end -mutable struct ModelInternal - container::OptimizationContainer - ic_model_container::Union{Nothing, OptimizationContainer} +mutable struct ModelInternal{T <: AbstractModelContainer} + container::T + ic_model_container::Union{Nothing, T} status::BuildStatus run_status::RunStatus base_conversion::Bool @@ -31,11 +31,11 @@ mutable struct ModelInternal end function ModelInternal( - container::OptimizationContainer; + container::T; ext = Dict{String, Any}(), recorders = [], -) - return ModelInternal( +) where {T <: AbstractModelContainer} + return ModelInternal{T}( container, nothing, BuildStatus.EMPTY, diff --git a/src/operation/operation_model_interface.jl b/src/operation/operation_model_interface.jl index 1677dbfab3..faa323ae6b 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -10,7 +10,7 @@ get_execution_count(model::OperationModel) = get_internal(model).execution_count get_executions(model::OperationModel) = get_internal(model).executions get_initial_time(model::OperationModel) = get_initial_time(get_settings(model)) get_internal(model::OperationModel) = model.internal -get_jump_model(model::OperationModel) = get_internal(model).container.JuMPmodel +get_jump_model(model::OperationModel) = get_jump_model(get_internal(model).container) get_name(model::OperationModel) = model.name get_store(model::OperationModel) = model.store is_synchronized(model::OperationModel) = is_synchronized(get_optimization_container(model)) diff --git a/test/performance/performance_test.jl b/test/performance/performance_test.jl index f7b0145b29..9d62bc281d 100644 --- a/test/performance/performance_test.jl +++ b/test/performance/performance_test.jl @@ -23,7 +23,7 @@ try for i in 1:2 template_uc = ProblemTemplate( NetworkModel( - StandardPTDFModel; + PTDFPowerModel; use_slacks = true, PTDF_matrix = PTDF(sys_rts_da), duals = [CopperPlateBalanceConstraint], diff --git a/test/run_partitioned_simulation.jl b/test/run_partitioned_simulation.jl index 6cbb5684cb..75e2a99d2c 100644 --- a/test/run_partitioned_simulation.jl +++ b/test/run_partitioned_simulation.jl @@ -113,7 +113,7 @@ function build_simulation( set_network_model!( template_uc, NetworkModel( - StandardPTDFModel; + PTDFPowerModel; PTDF_matrix = PTDF(c_sys5_pjm_da), # duals = [CopperPlateBalanceConstraint] ), diff --git a/test/test_device_branch_constructors.jl b/test/test_device_branch_constructors.jl index c6d404ca80..6768d7dc1d 100644 --- a/test/test_device_branch_constructors.jl +++ b/test/test_device_branch_constructors.jl @@ -1,7 +1,7 @@ @testset "DC Power Flow Models Monitored Line Flow Constraints and Static Unbounded" begin system = PSB.build_system(PSITestSystems, "c_sys5_ml") limits = PSY.get_flow_limits(PSY.get_component(MonitoredLine, system, "1")) - for model in [DCPPowerModel, StandardPTDFModel] + for model in [DCPPowerModel, PTDFPowerModel] template = get_thermal_dispatch_template_network( NetworkModel(model; PTDF_matrix = PTDF(system)), ) @@ -47,7 +47,7 @@ end @testset "DC Power Flow Models Monitored Line Flow Constraints and Static with inequalities" begin system = PSB.build_system(PSITestSystems, "c_sys5_ml") set_rate!(PSY.get_component(Line, system, "2"), 1.5) - for model in [DCPPowerModel, StandardPTDFModel] + for model in [DCPPowerModel, PTDFPowerModel] template = get_thermal_dispatch_template_network( NetworkModel(model; PTDF_matrix = PTDF(system)), ) @@ -67,7 +67,7 @@ end @testset "DC Power Flow Models Monitored Line Flow Constraints and Static with Bounds" begin system = PSB.build_system(PSITestSystems, "c_sys5_ml") set_rate!(PSY.get_component(Line, system, "2"), 1.5) - for model in [DCPPowerModel, StandardPTDFModel] + for model in [DCPPowerModel, PTDFPowerModel] template = get_thermal_dispatch_template_network( NetworkModel(model; PTDF_matrix = PTDF(system)), ) @@ -106,7 +106,7 @@ end transformer = PSY.get_component(Transformer2W, system, "Trans4") rate_limit2w = PSY.get_rate(tap_transformer) - for model in [DCPPowerModel, StandardPTDFModel] + for model in [DCPPowerModel, PTDFPowerModel] template = get_template_dispatch_with_network( NetworkModel(model; PTDF_matrix = PTDF(system)), ) @@ -163,7 +163,7 @@ end transformer = PSY.get_component(Transformer2W, system, "Trans4") rate_limit2w = PSY.get_rate(tap_transformer) - for model in [DCPPowerModel, StandardPTDFModel] + for model in [DCPPowerModel, PTDFPowerModel] template = get_template_dispatch_with_network( NetworkModel(model; PTDF_matrix = PTDF(system)), ) @@ -237,7 +237,7 @@ end add_component!(sys_5, hvdc) template_uc = ProblemTemplate( - NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF(sys_5)), + NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(sys_5)), ) set_device_model!(template_uc, ThermalStandard, ThermalCompactUnitCommitment) @@ -311,7 +311,7 @@ end ) add_component!(sys_5, hvdc) - for net_model in [DCPPowerModel, StandardPTDFModel] + for net_model in [DCPPowerModel, PTDFPowerModel] @testset "$net_model" begin PSY.set_loss!(hvdc, (l0 = 0.0, l1 = 0.0)) template_uc = ProblemTemplate( @@ -495,7 +495,7 @@ end rate_limit2w = PSY.get_rate(tap_transformer) template = get_template_dispatch_with_network( - NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF(system)), + NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(system)), ) set_device_model!(template, DeviceModel(TwoTerminalHVDCLine, HVDCTwoTerminalLossless)) model_m = DecisionModel(template, system; optimizer = HiGHS_optimizer) @@ -555,7 +555,7 @@ end add_component!(system, ps) template = get_template_dispatch_with_network( - NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF(system)), + NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(system)), ) set_device_model!(template, DeviceModel(PhaseShiftingTransformer, PhaseAngleControl)) model_m = DecisionModel(template, system; optimizer = HiGHS_optimizer) diff --git a/test/test_device_hvdc.jl b/test/test_device_hvdc.jl index 42ce383876..cbc2980b4c 100644 --- a/test/test_device_hvdc.jl +++ b/test/test_device_hvdc.jl @@ -15,11 +15,11 @@ set_device_model!(template_uc, DeviceModel(TModelHVDCLine, LossLessLine)) model = DecisionModel(template_uc, sys_5; name = "UC", optimizer = HiGHS_optimizer) @test build!(model; output_dir = mktempdir()) == PSI.BuildStatus.BUILT - moi_tests(model, 1656, 0, 1536, 816, 888, true) + moi_tests(model, 1656, 288, 1248, 528, 888, true) @test solve!(model) == RunStatus.SUCCESSFUL template_uc = ProblemTemplate(NetworkModel( - StandardPTDFModel; + PTDFPowerModel; #use_slacks=true, PTDF_matrix = PTDF(sys_5), #duals=[CopperPlateBalanceConstraint], diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index ae77e6e3a1..4f9f1f931d 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -643,7 +643,7 @@ end c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") c_sys5_dc = PSB.build_system(PSITestSystems, "c_sys5_dc") systems = [c_sys5, c_sys5_dc] - networks = [DCPPowerModel, NFAPowerModel, StandardPTDFModel, CopperPlatePowerModel] + networks = [DCPPowerModel, NFAPowerModel, PTDFPowerModel, CopperPlatePowerModel] commitment_models = [ThermalStandardUnitCommitment, ThermalCompactUnitCommitment] PTDF_ref = IdDict{System, PTDF}(c_sys5 => PTDF(c_sys5), c_sys5_dc => PTDF(c_sys5_dc)) @@ -793,7 +793,7 @@ end @testset "Test Must Run ThermalGen" begin sys_5 = build_system(PSITestSystems, "c_sys5_uc") template_uc = - ProblemTemplate(NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF(sys_5))) + ProblemTemplate(NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(sys_5))) set_device_model!(template_uc, ThermalStandard, ThermalCompactUnitCommitment) set_device_model!(template_uc, RenewableDispatch, FixedOutput) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) diff --git a/test/test_model_decision.jl b/test/test_model_decision.jl index 13fe6803f4..5b9b91a7f8 100644 --- a/test/test_model_decision.jl +++ b/test/test_model_decision.jl @@ -122,7 +122,7 @@ end @testset "Decision Model Solve with Slacks" begin c_sys5_re = PSB.build_system(PSITestSystems, "c_sys5_re") - networks = [StandardPTDFModel, DCPPowerModel, ACPPowerModel] + networks = [PTDFPowerModel, DCPPowerModel, ACPPowerModel] for network in networks template = get_thermal_dispatch_template_network( NetworkModel(network; use_slacks = true, PTDF_matrix = PTDF(c_sys5_re)), @@ -154,8 +154,8 @@ end @test UC_output == RunStatus.SUCCESSFUL end -@testset "Test Locational Marginal Prices between DC lossless with PowerModels vs StandardPTDFModel" begin - networks = [DCPPowerModel, StandardPTDFModel] +@testset "Test Locational Marginal Prices between DC lossless with PowerModels vs PTDFPowerModel" begin + networks = [DCPPowerModel, PTDFPowerModel] sys = PSB.build_system(PSITestSystems, "c_sys5") ptdf = PTDF(sys) # These are the duals of interest for the test @@ -165,7 +165,7 @@ end template = get_template_dispatch_with_network( NetworkModel(network; PTDF_matrix = ptdf, duals = dual_constraint[ix]), ) - if network == StandardPTDFModel + if network == PTDFPowerModel set_device_model!( template, DeviceModel(PSY.Line, PSI.StaticBranch; duals = [NetworkFlowConstraint]), @@ -178,7 +178,7 @@ end res = ProblemResults(model) # These tests require results to be working - if network == StandardPTDFModel + if network == PTDFPowerModel push!(LMPs, abs.(psi_ptdf_lmps(res, ptdf))) else duals = read_dual(res, NodalBalanceActiveConstraint, ACBus) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 3f44b92848..817fe1381f 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -1,31 +1,30 @@ # Note to devs. Use GLPK or Cbc for models with linear constraints and linear cost functions # Use OSQP for models with quadratic cost function and linear constraints and ipopt otherwise -const networks_for_testing = - networks = [ - (PM.ACPPowerModel, fast_ipopt_optimizer), - (PM.ACRPowerModel, fast_ipopt_optimizer), - (PM.ACTPowerModel, fast_ipopt_optimizer), - #(PM.IVRPowerModel, fast_ipopt_optimizer), #instantiate_ivp_expr_model not implemented - (PM.DCPPowerModel, fast_ipopt_optimizer), - (PM.DCMPPowerModel, fast_ipopt_optimizer), - (PM.NFAPowerModel, fast_ipopt_optimizer), - (PM.DCPLLPowerModel, fast_ipopt_optimizer), - (PM.LPACCPowerModel, fast_ipopt_optimizer), - (PM.SOCWRPowerModel, fast_ipopt_optimizer), - (PM.SOCWRConicPowerModel, scs_solver), - (PM.QCRMPowerModel, fast_ipopt_optimizer), - (PM.QCLSPowerModel, fast_ipopt_optimizer), - #(PM.SOCBFPowerModel, fast_ipopt_optimizer), # not implemented - (PM.BFAPowerModel, fast_ipopt_optimizer), - #(PM.SOCBFConicPowerModel, fast_ipopt_optimizer), # not implemented - (PM.SDPWRMPowerModel, scs_solver), - (PM.SparseSDPWRMPowerModel, scs_solver), - (PTDFPowerModel, fast_ipopt_optimizer), - ] +const NETWORKS_FOR_TESTING = [ + (PM.ACPPowerModel, fast_ipopt_optimizer), + (PM.ACRPowerModel, fast_ipopt_optimizer), + (PM.ACTPowerModel, fast_ipopt_optimizer), + #(PM.IVRPowerModel, fast_ipopt_optimizer), #instantiate_ivp_expr_model not implemented + (PM.DCPPowerModel, fast_ipopt_optimizer), + (PM.DCMPPowerModel, fast_ipopt_optimizer), + (PM.NFAPowerModel, fast_ipopt_optimizer), + (PM.DCPLLPowerModel, fast_ipopt_optimizer), + (PM.LPACCPowerModel, fast_ipopt_optimizer), + (PM.SOCWRPowerModel, fast_ipopt_optimizer), + (PM.SOCWRConicPowerModel, scs_solver), + (PM.QCRMPowerModel, fast_ipopt_optimizer), + (PM.QCLSPowerModel, fast_ipopt_optimizer), + #(PM.SOCBFPowerModel, fast_ipopt_optimizer), # not implemented + (PM.BFAPowerModel, fast_ipopt_optimizer), + #(PM.SOCBFConicPowerModel, fast_ipopt_optimizer), # not implemented + (PM.SDPWRMPowerModel, scs_solver), + (PM.SparseSDPWRMPowerModel, scs_solver), +] + @testset "All PowerModels models construction" begin c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") - for (network, solver) in networks_for_testing + for (network, solver) in NETWORKS_FOR_TESTING template = get_thermal_dispatch_template_network( NetworkModel(network; PTDF_matrix = PTDF(c_sys5)), ) @@ -89,7 +88,7 @@ end end @testset "Network DC-PF with PTDF Model" begin - template = get_thermal_dispatch_template_network(StandardPTDFModel) + template = get_thermal_dispatch_template_network(PTDFPowerModel) c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") c_sys14 = PSB.build_system(PSITestSystems, "c_sys14") c_sys14_dc = PSB.build_system(PSITestSystems, "c_sys14_dc") @@ -118,7 +117,7 @@ end ) for (ix, sys) in enumerate(systems) template = get_thermal_dispatch_template_network( - NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF_ref[sys]), + NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF_ref[sys]), ) ps_model = DecisionModel(template, sys; optimizer = HiGHS_optimizer) @@ -179,58 +178,10 @@ end c_sys14_dc => 142000.0, ) for (ix, sys) in enumerate(systems) - template = get_thermal_dispatch_template_network(StandardPTDFModel) + template = get_thermal_dispatch_template_network(PTDFPowerModel) template = get_thermal_dispatch_template_network( - NetworkModel(StandardPTDFModel; PTDF_matrix = PTDF_ref[sys]), - ) - ps_model = DecisionModel(template, sys; optimizer = HiGHS_optimizer) - - @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT - psi_constraint_test(ps_model, constraint_keys) - moi_tests( - ps_model, - test_results[sys][1], - test_results[sys][2], - test_results[sys][3], - test_results[sys][4], - test_results[sys][5], - false, - ) - psi_checkobjfun_test(ps_model, objfuncs[ix]) - psi_checksolve_test( - ps_model, - [MOI.OPTIMAL, MOI.ALMOST_OPTIMAL], - test_obj_values[sys], - 10000, + NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF_ref[sys]), ) - end -end - -@testset "Sparse Network DC-PF with PTDFPowerModel" begin - c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") - c_sys14 = PSB.build_system(PSITestSystems, "c_sys14") - c_sys14_dc = PSB.build_system(PSITestSystems, "c_sys14_dc") - systems = [c_sys5, c_sys14, c_sys14_dc] - objfuncs = [GAEVF, GQEVF, GQEVF] - constraint_keys = [ - PSI.ConstraintKey(RateLimitConstraint, PSY.Line, "lb"), - PSI.ConstraintKey(RateLimitConstraint, PSY.Line, "ub"), - PSI.ConstraintKey(CopperPlateBalanceConstraint, PSY.System), - PSI.ConstraintKey(NetworkFlowConstraint, PSY.Line), - ] - test_results = IdDict{System, Vector{Int}}( - c_sys5 => [264, 0, 264, 264, 168], - c_sys14 => [600, 0, 600, 600, 504], - c_sys14_dc => [600, 0, 648, 552, 456], - ) - test_obj_values = IdDict{System, Float64}( - c_sys5 => 340000.0, - c_sys14 => 142000.0, - c_sys14_dc => 142000.0, - ) - for (ix, sys) in enumerate(systems) - template = get_thermal_dispatch_template_network(PTDFPowerModel) ps_model = DecisionModel(template, sys; optimizer = HiGHS_optimizer) @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == @@ -267,9 +218,9 @@ end PSI.ConstraintKey(PSI.NodalBalanceActiveConstraint, PSY.ACBus), ] test_results = IdDict{System, Vector{Int}}( - c_sys5 => [384, 0, 408, 408, 288], - c_sys14 => [936, 0, 1080, 1080, 840], - c_sys14_dc => [984, 0, 1080, 984, 840], + c_sys5 => [384, 144, 264, 264, 288], + c_sys14 => [936, 480, 600, 600, 840], + c_sys14_dc => [984, 432, 648, 552, 840], ) test_obj_values = IdDict{System, Float64}( c_sys5 => 342000.0, @@ -315,9 +266,9 @@ end PSI.ConstraintKey(PSI.NodalBalanceReactiveConstraint, PSY.ACBus), ] test_results = IdDict{System, Vector{Int}}( - c_sys5 => [1056, 0, 384, 384, 264], - c_sys14 => [2832, 0, 720, 720, 696], - c_sys14_dc => [2832, 0, 768, 672, 744], + c_sys5 => [1056, 144, 240, 240, 264], + c_sys14 => [2832, 480, 240, 240, 696], + c_sys14_dc => [2832, 432, 336, 240, 744], ) test_obj_values = IdDict{System, Float64}( c_sys5 => 340000.0, @@ -411,9 +362,9 @@ end c_sys14_dc => [2832, 0, 336, 240, 744], ) ACT_test_results = Dict{System, Vector{Int}}( - c_sys5 => [1344, 0, 384, 384, 840], - c_sys14 => [3792, 0, 720, 720, 2616], - c_sys14_dc => [3696, 0, 768, 672, 2472], + c_sys5 => [1344, 144, 240, 240, 840], + c_sys14 => [3792, 480, 240, 240, 2616], + c_sys14_dc => [3696, 432, 336, 240, 2472], ) test_results = Dict(zip(networks, [ACR_test_results, ACT_test_results])) for network in networks, sys in systems @@ -450,14 +401,14 @@ end c_sys14_dc => 142000.0, ) DCPLL_test_results = Dict{System, Vector{Int}}( - c_sys5 => [528, 0, 408, 408, 288], - c_sys14 => [1416, 0, 1080, 1080, 840], - c_sys14_dc => [1416, 0, 1080, 984, 840], + c_sys5 => [528, 144, 264, 264, 288], + c_sys14 => [1416, 480, 600, 600, 840], + c_sys14_dc => [1416, 432, 648, 552, 840], ) LPACC_test_results = Dict{System, Vector{Int}}( - c_sys5 => [1200, 0, 384, 384, 840], - c_sys14 => [3312, 0, 720, 720, 2616], - c_sys14_dc => [3264, 0, 768, 672, 2472], + c_sys5 => [1200, 144, 240, 240, 840], + c_sys14 => [3312, 480, 240, 240, 2616], + c_sys14_dc => [3264, 432, 336, 240, 2472], ) test_results = Dict(zip(networks, [DCPLL_test_results, LPACC_test_results])) for network in networks, (ix, sys) in enumerate(systems) @@ -563,7 +514,7 @@ end set_active_power_limits_to!(hvdc_link, (min = 0.0, max = 0.0)) # Test not passing the PTDF to the Template - template = get_thermal_dispatch_template_network(NetworkModel(StandardPTDFModel)) + template = get_thermal_dispatch_template_network(NetworkModel(PTDFPowerModel)) ps_model = DecisionModel(template, c_sys5; optimizer = HiGHS_optimizer) @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == PSI.BuildStatus.BUILT @@ -603,7 +554,7 @@ end c_sys5 = PSB.build_system(PSISystems, "2Area 5 Bus System") # Test passing a VirtualPTDF Model template = get_thermal_dispatch_template_network( - NetworkModel(StandardPTDFModel; PTDF_matrix = VirtualPTDF(c_sys5)), + NetworkModel(PTDFPowerModel; PTDF_matrix = VirtualPTDF(c_sys5)), ) ps_model = DecisionModel(template, c_sys5; optimizer = HiGHS_optimizer) @@ -663,7 +614,7 @@ end set_active_power_limits_to!(hvdc_link, (min = 0.0, max = 0.0)) # Test not passing the PTDF to the Template - template = get_thermal_dispatch_template_network(NetworkModel(StandardPTDFModel)) + template = get_thermal_dispatch_template_network(NetworkModel(PTDFPowerModel)) ps_model = DecisionModel(template, c_sys5; optimizer = HiGHS_optimizer) @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == PSI.BuildStatus.BUILT @@ -699,249 +650,82 @@ end @test all(isapprox.(sum(zone_2_gen .+ zone_2_load; dims = 2), 0.0; atol = 1e-3)) end -function _updated_5bus_sys_with_extensions() - sys = PSB.build_system(PSITestSystems, "c_sys5_uc") - new_sys = deepcopy(sys) - ################################ - #### Create Extension Buses #### - ################################ - - busC = get_component(ACBus, new_sys, "nodeC") - - busC_ext1 = ACBus(; - number = 301, - name = "nodeC_ext1", - bustype = ACBusTypes.PQ, - angle = 0.0, - magnitude = 1.0, - voltage_limits = (min = 0.9, max = 1.05), - base_voltage = 230.0, - area = nothing, - load_zone = nothing, - ) - - busC_ext2 = ACBus(; - number = 302, - name = "nodeC_ext2", - bustype = ACBusTypes.PQ, - angle = 0.0, - magnitude = 1.0, - voltage_limits = (min = 0.9, max = 1.05), - base_voltage = 230.0, - area = nothing, - load_zone = nothing, - ) - - add_components!(new_sys, [busC_ext1, busC_ext2]) - - ################################ - #### Create Extension Lines #### - ################################ - - line_C_to_ext1 = Line(; - name = "C_to_ext1", - available = true, - active_power_flow = 0.0, - reactive_power_flow = 0.0, - arc = Arc(; from = busC, to = busC_ext1), - #r = 0.00281, - r = 0.0, - x = 0.0281, - b = (from = 0.00356, to = 0.00356), - rate = 2.0, - angle_limits = (min = -0.7, max = 0.7), - ) - - line_ext1_to_ext2 = Line(; - name = "ext1_to_ext2", - available = true, - active_power_flow = 0.0, - reactive_power_flow = 0.0, - arc = Arc(; from = busC_ext1, to = busC_ext2), - #r = 0.00281, - r = 0.0, - x = 0.0281, - b = (from = 0.00356, to = 0.00356), - rate = 2.0, - angle_limits = (min = -0.7, max = 0.7), - ) - - add_components!(new_sys, [line_C_to_ext1, line_ext1_to_ext2]) - - ################################### - ###### Update Extension Loads ##### - ################################### - - load_bus3 = get_component(PowerLoad, new_sys, "Bus3") - - load_ext1 = PowerLoad(; - name = "Bus_ext1", - available = true, - bus = busC_ext1, - active_power = 1.0, - reactive_power = 0.9861 / 3, - base_power = 100.0, - max_active_power = 1.0, - max_reactive_power = 0.9861 / 3, - ) - - load_ext2 = PowerLoad(; - name = "Bus_ext2", - available = true, - bus = busC_ext2, - active_power = 1.0, - reactive_power = 0.9861 / 3, - base_power = 100.0, - max_active_power = 1.0, - max_reactive_power = 0.9861 / 3, - ) - - add_components!(new_sys, [load_ext1, load_ext2]) - - copy_time_series!(load_ext1, load_bus3) - copy_time_series!(load_ext2, load_bus3) - - set_active_power!(load_bus3, 1.0) - set_max_active_power!(load_bus3, 1.0) - set_reactive_power!(load_bus3, 0.3287) - set_max_reactive_power!(load_bus3, 0.3287) - return new_sys -end - -@testset "StandardPTDF Radial Branches Test" begin - new_sys = _updated_5bus_sys_with_extensions() - - net_model = StandardPTDFModel - - template_uc = template_unit_commitment(; - network = NetworkModel(net_model; - reduce_radial_branches = true, - use_slacks = false, - ), - ) - thermal_model = ThermalStandardUnitCommitment - set_device_model!(template_uc, ThermalStandard, thermal_model) - - ##### Solve Reduced Model #### - solver = GLPK_optimizer - uc_model_red = DecisionModel( - template_uc, - new_sys; - optimizer = solver, - name = "UC_RED", - store_variable_names = true, - ) - - @test build!(uc_model_red; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT - solve!(uc_model_red) - - res_red = ProblemResults(uc_model_red) - - flow_lines = read_variable(res_red, "FlowActivePowerVariable__Line") - line_names = DataFrames.names(flow_lines)[2:end] - - ##### Solve Original Model #### - template_uc_orig = template_unit_commitment(; - network = NetworkModel(net_model; - reduce_radial_branches = false, - use_slacks = false, - ), - ) - set_device_model!(template_uc_orig, ThermalStandard, thermal_model) - - uc_model_orig = DecisionModel( - template_uc_orig, - new_sys; - optimizer = solver, - name = "UC_ORIG", - store_variable_names = true, - ) - - @test build!(uc_model_orig; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT - solve!(uc_model_orig) - - res_orig = ProblemResults(uc_model_orig) - - flow_lines_orig = read_variable(res_orig, "FlowActivePowerVariable__Line") - - for line in line_names - @test isapprox(flow_lines[!, line], flow_lines_orig[!, line]) - end -end - -@testset "DCPPowerModel Radial Branches Test" begin - net_model = DCPPowerModel - - template_uc = template_unit_commitment(; - network = NetworkModel(net_model; - reduce_radial_branches = true, - use_slacks = false, - ), - ) - thermal_model = ThermalStandardUnitCommitment - set_device_model!(template_uc, ThermalStandard, thermal_model) - - ##### Solve Reduced Model #### - solver = GLPK_optimizer - uc_model_red = DecisionModel( - template_uc, - new_sys; - optimizer = solver, - name = "UC_RED", - store_variable_names = true, - ) - - @test build!(uc_model_red; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT - solve!(uc_model_red) +# These models are easier to test due to their lossless nature +@testset "StandardPTDF/DCPPowerModel Radial Branches Test" begin + new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") + for net_model in [DCPPowerModel, PTDFPowerModel] + template_uc = template_unit_commitment(; + network = NetworkModel(net_model; + reduce_radial_branches = true, + use_slacks = false, + ), + ) + thermal_model = ThermalStandardUnitCommitment + set_device_model!(template_uc, ThermalStandard, thermal_model) + + ##### Solve Reduced Model #### + solver = GLPK_optimizer + uc_model_red = DecisionModel( + template_uc, + new_sys; + optimizer = solver, + name = "UC_RED", + store_variable_names = true, + ) - res_red = ProblemResults(uc_model_red) + @test build!(uc_model_red; output_dir = mktempdir(; cleanup = true)) == + PSI.BuildStatus.BUILT + solve!(uc_model_red) - flow_lines = read_variable(res_red, "FlowActivePowerVariable__Line") - line_names = DataFrames.names(flow_lines)[2:end] + res_red = ProblemResults(uc_model_red) - ##### Solve Original Model #### - template_uc_orig = template_unit_commitment(; - network = NetworkModel(net_model; - reduce_radial_branches = false, - use_slacks = false, - ), - ) - set_device_model!(template_uc_orig, ThermalStandard, thermal_model) + flow_lines = read_variable(res_red, "FlowActivePowerVariable__Line") + line_names = DataFrames.names(flow_lines)[2:end] - uc_model_orig = DecisionModel( - template_uc_orig, - new_sys; - optimizer = solver, - name = "UC_ORIG", - store_variable_names = true, - ) + ##### Solve Original Model #### + template_uc_orig = template_unit_commitment(; + network = NetworkModel(net_model; + reduce_radial_branches = false, + use_slacks = false, + ), + ) + set_device_model!(template_uc_orig, ThermalStandard, thermal_model) + + uc_model_orig = DecisionModel( + template_uc_orig, + new_sys; + optimizer = solver, + name = "UC_ORIG", + store_variable_names = true, + ) - @test build!(uc_model_orig; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT - solve!(uc_model_orig) + @test build!(uc_model_orig; output_dir = mktempdir(; cleanup = true)) == + PSI.BuildStatus.BUILT + solve!(uc_model_orig) - res_orig = ProblemResults(uc_model_orig) + res_orig = ProblemResults(uc_model_orig) - flow_lines_orig = read_variable(res_orig, "FlowActivePowerVariable__Line") + flow_lines_orig = read_variable(res_orig, "FlowActivePowerVariable__Line") - for line in line_names - @test isapprox(flow_lines[!, line], flow_lines_orig[!, line]) + for line in line_names + @test isapprox(flow_lines[!, line], flow_lines_orig[!, line]) + end end end @testset "All PowerModels models construction with reduced radial branches" begin - c_sys5 = _updated_5bus_sys_with_extensions() - for (network, solver) in networks_for_testing + new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") + for (network, solver) in NETWORKS_FOR_TESTING + if network ∈ PSI.INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS + continue + end template = get_thermal_dispatch_template_network( NetworkModel(network; - PTDF_matrix = PTDF(c_sys5), + PTDF_matrix = PTDF(new_sys), reduce_radial_branches = true, use_slacks = true), ) - ps_model = DecisionModel(template, c_sys5; optimizer = solver) + ps_model = DecisionModel(template, new_sys; optimizer = solver) @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == PSI.BuildStatus.BUILT @test ps_model.internal.container.pm !== nothing diff --git a/test/test_utils/operations_problem_templates.jl b/test/test_utils/operations_problem_templates.jl index 43c1e5bd1a..3e0a400f29 100644 --- a/test/test_utils/operations_problem_templates.jl +++ b/test/test_utils/operations_problem_templates.jl @@ -63,7 +63,7 @@ function get_template_hydro_st_ed(network = CopperPlatePowerModel, duals = []) return template end -function get_template_dispatch_with_network(network = StandardPTDFModel) +function get_template_dispatch_with_network(network = PTDFPowerModel) template = ProblemTemplate(network) set_device_model!(template, PowerLoad, StaticPowerLoad) set_device_model!(template, ThermalStandard, ThermalBasicDispatch)