From 54d707fb700312e26c82681461bb1a65194b5bb6 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 16:03:51 -0700 Subject: [PATCH 01/59] add radial_lines to network models --- src/core/network_model.jl | 21 +++++++++++++++++++-- 1 file changed, 19 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 058d512f77..adef1a9a3f 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -32,16 +32,27 @@ 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 + remove_radial_lines::Bool function NetworkModel( ::Type{T}; use_slacks = false, PTDF_matrix = nothing, + remove_radial_lines = false, subnetworks = Dict{Int, Set{Int}}(), duals = Vector{DataType}(), ) where {T <: PM.AbstractPowerModel} _check_pm_formulation(T) - new{T}(use_slacks, PTDF_matrix, subnetworks, Dict{PSY.ACBus, Int}(), duals) + new{T}( + use_slacks, + PTDF_matrix, + subnetworks, + Dict{PSY.ACBus, Int}(), + duals, + PNM.RadialBranches(), + remove_radial_lines, + ) end end @@ -69,6 +80,9 @@ function instantiate_network_model( if isempty(model.subnetworks) model.subnetworks = PNM.find_subnetworks(sys) end + if model.remove_radial_lines + model.radial_branches = PNM.RadialBranches(sys) + end return end @@ -90,7 +104,7 @@ end function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys::PSY.System) if get_PTDF_matrix(model) === nothing @info "PTDF Matrix not provided. Calculating using PowerNetworkMatrices.PTDF" - model.PTDF_matrix = PNM.PTDF(sys) + model.PTDF_matrix = PNM.PTDF(sys; remove_radial_lines = model.remove_radial_lines) end get_PTDF_matrix(model).subnetworks model.subnetworks = deepcopy(get_PTDF_matrix(model).subnetworks) @@ -98,6 +112,9 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: @debug "System Contains Multiple Subnetworks. Assigning buses to subnetworks." _assign_subnetworks_to_buses(model, sys) end + if model.remove_radial_lines + model.radial_branches = model.PTDF_matrix.radial_branches + end return end From df8f36115b937a1af49554e4e884ecbe88fb7d8e Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 16:52:49 -0700 Subject: [PATCH 02/59] add radial branches to model --- src/core/network_model.jl | 13 +++++++------ src/core/optimization_container.jl | 3 +-- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index adef1a9a3f..32a6f2be32 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,13 +33,13 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches - remove_radial_lines::Bool + reduce_radial_branches::Bool function NetworkModel( ::Type{T}; use_slacks = false, PTDF_matrix = nothing, - remove_radial_lines = false, + reduce_radial_branches = false, subnetworks = Dict{Int, Set{Int}}(), duals = Vector{DataType}(), ) where {T <: PM.AbstractPowerModel} @@ -51,7 +51,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), - remove_radial_lines, + reduce_radial_branches, ) end end @@ -80,7 +80,7 @@ function instantiate_network_model( if isempty(model.subnetworks) model.subnetworks = PNM.find_subnetworks(sys) end - if model.remove_radial_lines + if model.reduce_radial_branches model.radial_branches = PNM.RadialBranches(sys) end return @@ -104,7 +104,8 @@ end function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys::PSY.System) if get_PTDF_matrix(model) === nothing @info "PTDF Matrix not provided. Calculating using PowerNetworkMatrices.PTDF" - model.PTDF_matrix = PNM.PTDF(sys; remove_radial_lines = model.remove_radial_lines) + model.PTDF_matrix = + PNM.PTDF(sys; reduce_radial_branches = model.reduce_radial_branches) end get_PTDF_matrix(model).subnetworks model.subnetworks = deepcopy(get_PTDF_matrix(model).subnetworks) @@ -112,7 +113,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: @debug "System Contains Multiple Subnetworks. Assigning buses to subnetworks." _assign_subnetworks_to_buses(model, sys) end - if model.remove_radial_lines + if model.reduce_radial_branches model.radial_branches = model.PTDF_matrix.radial_branches end return diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index be85ae1a9c..89ba392b50 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -519,8 +519,7 @@ function build_impl!( container, transmission, transmission_model.subnetworks, - sys, - ) + sys,) # Order is required for device_model in values(template.devices) From fdba9badae58680ed5283da7d335f30dbc25f6d2 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 17:41:01 -0700 Subject: [PATCH 03/59] add bus maps --- src/core/network_model.jl | 2 ++ src/core/optimization_container.jl | 28 +++++++++++++++++++++++----- 2 files changed, 25 insertions(+), 5 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 32a6f2be32..8c5a21abae 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,6 +33,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches + bus_aggregation_map::Dict{Int, Int} reduce_radial_branches::Bool function NetworkModel( @@ -51,6 +52,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), + Dict{Int, Int}(), reduce_radial_branches, ) end diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 89ba392b50..82e5fb88b2 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -431,9 +431,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractPowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) - ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.ACBus) => _make_container_array(ac_bus_numbers, time_steps), @@ -450,9 +455,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractActivePowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) - ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.ACBus) => _make_container_array(ac_bus_numbers, time_steps), @@ -467,6 +477,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, ::Vector{Int}, ::Type{CopperPlatePowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) subnetworks_ref_buses = collect(keys(subnetworks)) @@ -482,9 +493,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{T}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) where {T <: Union{PTDFPowerModel, StandardPTDFModel}} time_steps = get_time_steps(container) - ac_bus_numbers = sort!(collect(Iterators.flatten(values(subnetworks)))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end subnetworks = collect(keys(subnetworks)) container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.System) => @@ -502,9 +518,10 @@ function initialize_system_expressions!( ::Type{T}, subnetworks::Dict{Int, Set{Int}}, system::PSY.System, + bus_reduction_map::Dict{Int64, Set{Int64}} ) where {T <: PM.AbstractPowerModel} dc_bus_numbers = [PSY.get_number(b) for b in PSY.get_components(PSY.DCBus, system)] - _make_system_expressions!(container, subnetworks, dc_bus_numbers, T) + _make_system_expressions!(container, subnetworks, dc_bus_numbers, T, bus_reduction_map) return end @@ -519,7 +536,8 @@ function build_impl!( container, transmission, transmission_model.subnetworks, - sys,) + sys, + transmission_model.radial_branches.bus_reduction_map) # Order is required for device_model in values(template.devices) From 040b15037d21ee7611eb67ffe1cb3ad577f3aad6 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 08:34:03 -0700 Subject: [PATCH 04/59] add missing sort! call --- src/core/optimization_container.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 82e5fb88b2..f6074b8e64 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -508,7 +508,7 @@ function _make_system_expressions!( ExpressionKey(ActivePowerBalance, PSY.DCBus) => _make_container_array(dc_bus_numbers, time_steps), ExpressionKey(ActivePowerBalance, PSY.ACBus) => - _make_container_array(ac_bus_numbers, time_steps), + _make_container_array(sort!(ac_bus_numbers), time_steps), ) return end From b8f96296374a29fdfe75a16e8ee1b3ce9dac573b Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:10 -0700 Subject: [PATCH 05/59] updates to network model --- src/core/network_model.jl | 5 +++-- src/initial_conditions/initialization.jl | 4 ++++ 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 8c5a21abae..15c0edbee6 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,7 +33,6 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches - bus_aggregation_map::Dict{Int, Int} reduce_radial_branches::Bool function NetworkModel( @@ -52,7 +51,6 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), - Dict{Int, Int}(), reduce_radial_branches, ) end @@ -60,6 +58,8 @@ 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_duals(m::NetworkModel) = m.duals get_network_formulation(::NetworkModel{T}) where {T} = T get_reference_buses(m::NetworkModel{T}) where {T <: PM.AbstractPowerModel} = @@ -116,6 +116,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: _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 diff --git a/src/initial_conditions/initialization.jl b/src/initial_conditions/initialization.jl index 80ff0f1fd7..d866993118 100644 --- a/src/initial_conditions/initialization.jl +++ b/src/initial_conditions/initialization.jl @@ -5,7 +5,11 @@ function get_initial_conditions_template(model::OperationModel) get_network_formulation(model.template); use_slacks = get_use_slacks(get_network_model(model.template)), PTDF_matrix = get_PTDF_matrix(get_network_model(model.template)), + reduce_radial_branches = get_reduce_radial_branches( + get_network_model(model.template), + ), ) + network_model.radial_branches = get_radial_branches(get_network_model(model.template)) network_model.subnetworks = get_subnetworks(get_network_model(model.template)) bus_area_map = get_bus_area_map(get_network_model(model.template)) if !isempty(bus_area_map) From 2deda19fb83eb3e1c3219641a5bfb87b373f9279 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:25 -0700 Subject: [PATCH 06/59] updates to the container --- src/core/optimization_container.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index f6074b8e64..5e68f6e067 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -431,7 +431,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractPowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -455,7 +455,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractActivePowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -477,7 +477,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, ::Vector{Int}, ::Type{CopperPlatePowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) subnetworks_ref_buses = collect(keys(subnetworks)) @@ -493,7 +493,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{T}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) where {T <: Union{PTDFPowerModel, StandardPTDFModel}} time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -518,7 +518,7 @@ function initialize_system_expressions!( ::Type{T}, subnetworks::Dict{Int, Set{Int}}, system::PSY.System, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) where {T <: PM.AbstractPowerModel} dc_bus_numbers = [PSY.get_number(b) for b in PSY.get_components(PSY.DCBus, system)] _make_system_expressions!(container, subnetworks, dc_bus_numbers, T, bus_reduction_map) From d18aaea8d5bc5d18fed4b9c8edb811721cc08557 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:33 -0700 Subject: [PATCH 07/59] updates to branch modeling --- src/devices_models/devices/AC_branches.jl | 49 ++++++++- .../devices/common/add_to_expression.jl | 99 +++++++++++++------ 2 files changed, 112 insertions(+), 36 deletions(-) diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 1176f83700..bf24ae3c3c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -168,14 +168,55 @@ function add_constraints!( container::OptimizationContainer, cons_type::Type{RateLimitConstraint}, devices::IS.FlattenIteratorWrapper{T}, - model::DeviceModel{T, U}, - ::NetworkModel{X}, + ::DeviceModel{T, U}, + network_model::NetworkModel{V}, ) where { T <: PSY.ACBranch, U <: AbstractBranchFormulation, - X <: PM.AbstractActivePowerModel, + V <: PM.AbstractActivePowerModel, } - add_range_constraints!(container, cons_type, FlowActivePowerVariable, devices, model, X) + time_steps = get_time_steps(container) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + device_names = [PSY.get_name(d) for d in devices] + else + device_names = PNM.get_meshed_branches(radial_branches) + end + + con_lb = + add_constraints_container!( + container, + cons_type(), + T, + device_names, + time_steps; + meta = "lb", + ) + con_ub = + add_constraints_container!( + container, + cons_type(), + T, + device_names, + time_steps; + meta = "ub", + ) + + array = get_variable(container, FlowActivePowerVariable(), T) + + for device in devices + ci_name = PSY.get_name(device) + if ci_name ∈ PNM.get_radial_branches(radial_branches) + continue + end + limits = get_min_max_limits(device, RateLimitConstraint, U) # depends on constraint type and formulation type + for t in time_steps + con_ub[ci_name, t] = + JuMP.@constraint(container.JuMPmodel, array[ci_name, t] <= limits.max) + con_lb[ci_name, t] = + JuMP.@constraint(container.JuMPmodel, array[ci_name, t] >= limits.min) + end + end return end diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 54aa14a590..22b56059b2 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -79,7 +79,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, model::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: TimeSeriesParameter, @@ -89,11 +89,12 @@ function add_to_expression!( } param_container = get_parameter(container, U(), V) multiplier = get_multiplier_array(param_container) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) name = PSY.get_name(d) _add_to_jump_expression!( - get_expression(container, T(), PSY.ACBus)[bus_number, t], + get_expression(container, T(), PSY.ACBus)[bus_no, t], get_parameter_column_refs(param_container, name)[t], multiplier[name, t], ) @@ -107,7 +108,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: ActivePowerBalance, U <: OnStatusParameter, @@ -116,13 +117,13 @@ function add_to_expression!( X <: PM.AbstractPowerModel, } parameter = get_parameter_array(container, U(), V) - + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) name = PSY.get_name(d) mult = get_expression_multiplier(U(), T(), d, W()) _add_to_jump_expression!( - get_expression(container, T(), PSY.ACBus)[bus_number, t], + get_expression(container, T(), PSY.ACBus)[bus_no, t], parameter[name, t], mult, ) @@ -139,7 +140,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: VariableType, @@ -149,11 +150,12 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -215,8 +217,10 @@ 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) for d in devices - bus_no_to = PSY.get_number(PSY.get_arc(d).to) + bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) + bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_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) @@ -250,8 +254,10 @@ 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) for d in devices - bus_no_from = PSY.get_number(PSY.get_arc(d).from) + bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) + bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_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) @@ -352,12 +358,14 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_arc(d).from) + bus_no_ = PSY.get_number(PSY.get_arc(d).from) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -385,12 +393,14 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_arc(d).to) + bus_no_ = PSY.get_number(PSY.get_arc(d).to) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -415,11 +425,13 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no_ = PSY.get_number(PSY.get_bus(d)) + PNM.get_mapped_bus_number(radial_branches, bus_no_) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), d, W()), ) @@ -574,10 +586,12 @@ function add_to_expression!( 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) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) param = get_parameter_column_refs(param_container, name) for t in get_time_steps(container) @@ -605,9 +619,11 @@ function add_to_expression!( 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) 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_ = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, 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) @@ -637,10 +653,13 @@ function add_to_expression!( 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) + @assert !isempty(radial_branches) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -675,10 +694,12 @@ function add_to_expression!( 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) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -715,16 +736,21 @@ function add_to_expression!( } var = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(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_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).from), t], + expression[bus_no_from, t], flow_variable, -1.0, ) _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).to), t], + expression[bus_no_to, t], flow_variable, 1.0, ) @@ -752,9 +778,12 @@ 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) for d in devices - bus_no_from = PSY.get_number(PSY.get_arc(d).from) - bus_no_to = PSY.get_number(PSY.get_arc(d).to) + 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_) 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) @@ -826,16 +855,21 @@ function add_to_expression!( ) 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) 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_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).from), t], + expression[bus_no_from, t], flow_variable, -get_variable_multiplier(U(), d, V()), ) _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).to), t], + expression[bus_no_to, t], flow_variable, get_variable_multiplier(U(), d, V()), ) @@ -1079,6 +1113,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus) expression = get_expression(container, T(), PSY.ACBus) @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(variable, 1) _add_to_jump_expression!( expression[n, t], @@ -1102,8 +1137,8 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus, "P") expression = get_expression(container, T(), PSY.ACBus) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) - for t in get_time_steps(container), n in bus_numbers + # 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(variable, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1126,8 +1161,8 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus, "Q") expression = get_expression(container, T(), PSY.ACBus) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) - for t in get_time_steps(container), n in bus_numbers + # 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(variable, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], From c393f8154de562d0420145971ef89999032f589a Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 16:32:53 -0700 Subject: [PATCH 08/59] remove debugging assertion --- src/devices_models/devices/common/add_to_expression.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 22b56059b2..f37e4e729c 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -654,7 +654,6 @@ function add_to_expression!( sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) radial_branches = get_radial_branches(network_model) - @assert !isempty(radial_branches) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) From c8845f159486dae58cd20b383b7316333b76e5bd Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 17:11:16 -0700 Subject: [PATCH 09/59] fix tests --- src/devices_models/devices/common/add_to_expression.jl | 6 +++--- test/test_utils/mock_operation_models.jl | 1 + 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index f37e4e729c..31340271dd 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -383,7 +383,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: ActivePowerBalance, U <: FlowActivePowerToFromVariable, @@ -415,7 +415,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: OnVariable, @@ -429,7 +429,7 @@ function add_to_expression!( for d in devices, t in get_time_steps(container) name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) - PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) _add_to_jump_expression!( expression[bus_no, t], variable[name, t], diff --git a/test/test_utils/mock_operation_models.jl b/test/test_utils/mock_operation_models.jl index a41128efb5..759684109c 100644 --- a/test/test_utils/mock_operation_models.jl +++ b/test/test_utils/mock_operation_models.jl @@ -117,6 +117,7 @@ function mock_construct_device!( PSI.get_network_formulation(template), PSI.get_network_model(template).subnetworks, PSI.get_system(problem), + Dict{Int64, Set{Int64}}(), ) if PSI.validate_available_devices(model, PSI.get_system(problem)) PSI.construct_device!( From 4d2c34cc33c720e58082a93f1f6fdbd7cfe08e42 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 21:26:19 -0700 Subject: [PATCH 10/59] fix bug with serialize function --- src/operation/operation_model_interface.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/operation/operation_model_interface.jl b/src/operation/operation_model_interface.jl index 688f6013f7..1677dbfab3 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -415,7 +415,7 @@ end function serialize_optimization_model(model::OperationModel, save_path::String) serialize_jump_optimization_model( - get_optimization_container(model), + get_jump_model(get_optimization_container(model)), save_path, ) return From 21b0ae40ff2fadd2d25a64c8ac0908cefd39097f Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 13:40:20 -0800 Subject: [PATCH 11/59] update branch constructors to ignore radial branches --- .../device_constructors/branch_constructor.jl | 8 ++-- src/devices_models/devices/AC_branches.jl | 37 ++++++++++++++++--- src/network_models/pm_translator.jl | 16 +++++--- src/network_models/powermodels_interface.jl | 35 +++++++++++++----- 4 files changed, 73 insertions(+), 23 deletions(-) diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index ec1b7d6d9d..616500cafb 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -106,7 +106,7 @@ function construct_device!( container, devices, device_model, - get_network_formulation(network_model), + network_model, ) end return @@ -273,7 +273,7 @@ function construct_device!( container, devices, device_model, - get_network_formulation(network_model), + network_model, ) add_constraint_dual!(container, sys, device_model) return @@ -331,7 +331,7 @@ function construct_device!( ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) - branch_rate_bounds!(container, devices, model, get_network_formulation(network_model)) + branch_rate_bounds!(container, devices, model, network_model) add_constraints!(container, RateLimitConstraintFromTo, devices, model, network_model) add_constraints!(container, RateLimitConstraintToFrom, devices, model, network_model) @@ -356,7 +356,7 @@ function construct_device!( ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) - branch_rate_bounds!(container, devices, model, get_network_formulation(network_model)) + branch_rate_bounds!(container, devices, model, network_model) add_constraint_dual!(container, sys, model) return end diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index bf24ae3c3c..de74ab915c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -104,11 +104,18 @@ function branch_rate_bounds!( container::OptimizationContainer, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::Type{<:PM.AbstractDCPModel}, + network_model::NetworkModel{<:PM.AbstractDCPModel}, ) 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) + for d in devices name = PSY.get_name(d) + if name ∈ radial_branches_names + continue + end for t in get_time_steps(container) JuMP.set_upper_bound(var[name, t], PSY.get_rate(d)) JuMP.set_lower_bound(var[name, t], -1.0 * PSY.get_rate(d)) @@ -121,15 +128,23 @@ function branch_rate_bounds!( container::OptimizationContainer, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::Type{<:PM.AbstractPowerModel}, + network_model::NetworkModel{<:PM.AbstractPowerModel}, ) where {B <: PSY.ACBranch} vars = [ get_variable(container, FlowActivePowerFromToVariable(), B), get_variable(container, FlowActivePowerToFromVariable(), B), ] + + time_steps = get_time_steps(container) + radial_branches = get_radial_branches(network_model) + radial_branches_names = PNM.get_radial_branches(radial_branches) + for d in devices name = PSY.get_name(d) - for t in get_time_steps(container), var in vars + if name ∈ radial_branches_names + continue + end + for t in time_steps, var in vars JuMP.set_upper_bound(var[name, t], PSY.get_rate(d)) JuMP.set_lower_bound(var[name, t], -1.0 * PSY.get_rate(d)) end @@ -256,7 +271,7 @@ function add_constraints!( cons_type::Type{RateLimitConstraintFromTo}, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::NetworkModel{T}, + network_model::NetworkModel{T}, ) where {B <: PSY.ACBranch, T <: PM.AbstractPowerModel} rating_data = [(PSY.get_name(h), PSY.get_rate(h)) for h in devices] @@ -272,7 +287,13 @@ 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) + for r in rating_data + if r[1] ∈ radial_branches_names + continue + end for t in time_steps constraint[r[1], t] = JuMP.@constraint( get_jump_model(container), @@ -291,7 +312,7 @@ function add_constraints!( cons_type::Type{RateLimitConstraintToFrom}, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::NetworkModel{T}, + network_model::NetworkModel{T}, ) where {B <: PSY.ACBranch, T <: PM.AbstractPowerModel} rating_data = [(PSY.get_name(h), PSY.get_rate(h)) for h in devices] @@ -307,7 +328,13 @@ 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) + for r in rating_data + if r[1] ∈ radial_branches_names + continue + end for t in time_steps constraint[r[1], t] = JuMP.@constraint( get_jump_model(container), diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 627a5b7cf0..ec012ff097 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -381,7 +381,7 @@ end function get_branches_to_pm( sys::PSY.System, - ::Type{S}, + network_model::NetworkModel{S}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -389,6 +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) for (d, device_model) in branch_template comp_type = get_component_type(device_model) if comp_type <: TwoTerminalHVDCTypes @@ -398,6 +400,10 @@ function get_branches_to_pm( start_idx += length(PM_branches) 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" + continue + end ix = i + start_idx PM_branches["$(ix)"] = get_branch_to_pm(ix, branch, get_formulation(device_model), S) @@ -413,7 +419,7 @@ end function get_branches_to_pm( sys::PSY.System, - ::Type{S}, + network_model::NetworkModel{S}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -442,7 +448,7 @@ end function get_branches_to_pm( ::PSY.System, - ::Type{PTDFPowerModel}, + network_model::NetworkModel{PTDFPowerModel}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -485,13 +491,13 @@ end function pass_to_pm(sys::PSY.System, template::ProblemTemplate, time_periods::Int) ac_lines, PMmap_ac = get_branches_to_pm( sys, - get_network_formulation(template), + get_network_model(template), PSY.ACBranch, template.branches, ) two_terminal_dc_lines, PMmap_dc = get_branches_to_pm( sys, - get_network_formulation(template), + get_network_model(template), TwoTerminalHVDCTypes, template.branches, length(ac_lines), diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index 932f2c158f..942a735985 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -291,17 +291,25 @@ function powermodels_network!( ) where {S <: PM.AbstractPowerModel} 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) - for t in time_steps, bus in buses - pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_p"] = + network_model = get_network_model(template) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_reduction_map = radial_branches.bus_reduction_map + ac_bus_numbers = collect(keys(bus_reduction_map)) + end + + for t in time_steps, bus_no in ac_bus_numbers + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_p"] = container.expressions[ExpressionKey(ActivePowerBalance, PSY.ACBus)][ - bus.number, + bus_no, t, ] - pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_q"] = + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_q"] = container.expressions[ExpressionKey(ReactivePowerBalance, PSY.ACBus)][ - bus.number, + bus_no, t, ] end @@ -324,10 +332,19 @@ function powermodels_network!( pm_data, PM_map = pass_to_pm(sys, template, time_steps[end]) buses = get_available_components(PSY.ACBus, sys) - for t in time_steps, bus in buses - pm_data["nw"]["$(t)"]["bus"]["$(PSY.get_number(bus))"]["inj_p"] = + network_model = get_network_model(template) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_reduction_map = radial_branches.bus_reduction_map + ac_bus_numbers = collect(keys(bus_reduction_map)) + end + + for t in time_steps, bus_no in ac_bus_numbers + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_p"] = container.expressions[ExpressionKey(ActivePowerBalance, PSY.ACBus)][ - PSY.get_number(bus), + bus_no, t, ] # pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_q"] = 0.0 From 8154c9afac7d0edc598caad98d91380e78112c94 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 13:40:40 -0800 Subject: [PATCH 12/59] update slack expressions --- src/devices_models/devices/common/add_to_expression.jl | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 31340271dd..17bdefdbb0 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -1111,9 +1111,10 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus) expression = get_expression(container, T(), PSY.ACBus) - @assert_op length(axes(variable, 1)) == length(axes(expression, 1)) + #TODO: Update Slacks to consider reduced buses + #@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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1137,7 +1138,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus, "P") expression = get_expression(container, T(), PSY.ACBus) # 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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1161,7 +1162,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus, "Q") expression = get_expression(container, T(), PSY.ACBus) # 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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], From d1dd6b1c1bdfba58e4561d2404a623fcfcfbe35a Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 15:28:08 -0800 Subject: [PATCH 13/59] add radial test --- test/test_radial_branches.jl | 234 +++++++++++++++++++++++++++++++++++ 1 file changed, 234 insertions(+) create mode 100644 test/test_radial_branches.jl diff --git a/test/test_radial_branches.jl b/test/test_radial_branches.jl new file mode 100644 index 0000000000..1ba2633443 --- /dev/null +++ b/test/test_radial_branches.jl @@ -0,0 +1,234 @@ +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 + new_sys = _updated_5bus_sys_with_extensions() + + 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) + + 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 From bf2407726bd6662c1261236a9202b5575e436ec4 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 8 Jan 2024 10:12:22 -0700 Subject: [PATCH 14/59] add new test --- test/test_network_constructors.jl | 62 ++++++++++++++++++++----------- 1 file changed, 40 insertions(+), 22 deletions(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 4de655700b..d6ffee3d1d 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -1,29 +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), +] + @testset "All PowerModels models construction" begin - 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), - ] c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") - for (network, solver) in networks + for (network, solver) in networks_for_testing template = get_thermal_dispatch_template_network( NetworkModel(network; PTDF_matrix = PTDF(c_sys5)), ) @@ -696,3 +697,20 @@ end ) @test all(isapprox.(sum(zone_2_gen .+ zone_2_load; dims = 2), 0.0; atol = 1e-3)) end + +@testset "All PowerModels models construction with reduced radial branches" begin + # use the system with a radial branch + #c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") + for (network, solver) in networks_for_testing + template = get_thermal_dispatch_template_network( + NetworkModel(network; + PTDF_matrix = PTDF(c_sys5), + reduce_radial_branches = true, + use_slacks), + ) + ps_model = DecisionModel(template, c_sys5; optimizer = solver) + @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == + PSI.BuildStatus.BUILT + @test ps_model.internal.container.pm !== nothing + end +end From f5e52bdc0b0e0da492b97588919ff10655c1de02 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 8 Jan 2024 10:13:24 -0700 Subject: [PATCH 15/59] update testing --- test/test_network_constructors.jl | 237 +++++++++++++++++++++++++++++- test/test_radial_branches.jl | 234 ----------------------------- 2 files changed, 235 insertions(+), 236 deletions(-) delete mode 100644 test/test_radial_branches.jl diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index d6ffee3d1d..0d30eb4456 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -698,9 +698,242 @@ 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) + + 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 "All PowerModels models construction with reduced radial branches" begin - # use the system with a radial branch - #c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") + c_sys5 = _updated_5bus_sys_with_extensions() for (network, solver) in networks_for_testing template = get_thermal_dispatch_template_network( NetworkModel(network; diff --git a/test/test_radial_branches.jl b/test/test_radial_branches.jl deleted file mode 100644 index 1ba2633443..0000000000 --- a/test/test_radial_branches.jl +++ /dev/null @@ -1,234 +0,0 @@ -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 - new_sys = _updated_5bus_sys_with_extensions() - - 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) - - 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 From 4f50cc8fe6fd10dd7c906e1f1da742772f1195ff Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 11 Jan 2024 10:15:55 -0700 Subject: [PATCH 16/59] bump PNM version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 92630cf604..223a9a623c 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,7 @@ LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" PowerModels = "~0.19" -PowerNetworkMatrices = "^0.9" +PowerNetworkMatrices = "^0.10" PowerSystems = "^3" PrettyTables = "2" ProgressMeter = "^1.5" From b50e02c60c5cd4eb230075342083909373d67787 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 11 Jan 2024 14:00:30 -0700 Subject: [PATCH 17/59] bump PNM dependency --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 223a9a623c..92630cf604 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,7 @@ LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" PowerModels = "~0.19" -PowerNetworkMatrices = "^0.10" +PowerNetworkMatrices = "^0.9" PowerSystems = "^3" PrettyTables = "2" ProgressMeter = "^1.5" From 5aefd0ad2d57e2c9aad907d9253ce6e73f1145d8 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 16:03:51 -0700 Subject: [PATCH 18/59] add radial_lines to network models --- src/core/network_model.jl | 21 +++++++++++++++++++-- 1 file changed, 19 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 058d512f77..adef1a9a3f 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -32,16 +32,27 @@ 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 + remove_radial_lines::Bool function NetworkModel( ::Type{T}; use_slacks = false, PTDF_matrix = nothing, + remove_radial_lines = false, subnetworks = Dict{Int, Set{Int}}(), duals = Vector{DataType}(), ) where {T <: PM.AbstractPowerModel} _check_pm_formulation(T) - new{T}(use_slacks, PTDF_matrix, subnetworks, Dict{PSY.ACBus, Int}(), duals) + new{T}( + use_slacks, + PTDF_matrix, + subnetworks, + Dict{PSY.ACBus, Int}(), + duals, + PNM.RadialBranches(), + remove_radial_lines, + ) end end @@ -69,6 +80,9 @@ function instantiate_network_model( if isempty(model.subnetworks) model.subnetworks = PNM.find_subnetworks(sys) end + if model.remove_radial_lines + model.radial_branches = PNM.RadialBranches(sys) + end return end @@ -90,7 +104,7 @@ end function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys::PSY.System) if get_PTDF_matrix(model) === nothing @info "PTDF Matrix not provided. Calculating using PowerNetworkMatrices.PTDF" - model.PTDF_matrix = PNM.PTDF(sys) + model.PTDF_matrix = PNM.PTDF(sys; remove_radial_lines = model.remove_radial_lines) end get_PTDF_matrix(model).subnetworks model.subnetworks = deepcopy(get_PTDF_matrix(model).subnetworks) @@ -98,6 +112,9 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: @debug "System Contains Multiple Subnetworks. Assigning buses to subnetworks." _assign_subnetworks_to_buses(model, sys) end + if model.remove_radial_lines + model.radial_branches = model.PTDF_matrix.radial_branches + end return end From 1ee92de9173ec711cfc2cd3ea82f9ffa012f5976 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 16:52:49 -0700 Subject: [PATCH 19/59] add radial branches to model --- src/core/network_model.jl | 13 +++++++------ src/core/optimization_container.jl | 3 +-- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index adef1a9a3f..32a6f2be32 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,13 +33,13 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches - remove_radial_lines::Bool + reduce_radial_branches::Bool function NetworkModel( ::Type{T}; use_slacks = false, PTDF_matrix = nothing, - remove_radial_lines = false, + reduce_radial_branches = false, subnetworks = Dict{Int, Set{Int}}(), duals = Vector{DataType}(), ) where {T <: PM.AbstractPowerModel} @@ -51,7 +51,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), - remove_radial_lines, + reduce_radial_branches, ) end end @@ -80,7 +80,7 @@ function instantiate_network_model( if isempty(model.subnetworks) model.subnetworks = PNM.find_subnetworks(sys) end - if model.remove_radial_lines + if model.reduce_radial_branches model.radial_branches = PNM.RadialBranches(sys) end return @@ -104,7 +104,8 @@ end function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys::PSY.System) if get_PTDF_matrix(model) === nothing @info "PTDF Matrix not provided. Calculating using PowerNetworkMatrices.PTDF" - model.PTDF_matrix = PNM.PTDF(sys; remove_radial_lines = model.remove_radial_lines) + model.PTDF_matrix = + PNM.PTDF(sys; reduce_radial_branches = model.reduce_radial_branches) end get_PTDF_matrix(model).subnetworks model.subnetworks = deepcopy(get_PTDF_matrix(model).subnetworks) @@ -112,7 +113,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: @debug "System Contains Multiple Subnetworks. Assigning buses to subnetworks." _assign_subnetworks_to_buses(model, sys) end - if model.remove_radial_lines + if model.reduce_radial_branches model.radial_branches = model.PTDF_matrix.radial_branches end return diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index be85ae1a9c..89ba392b50 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -519,8 +519,7 @@ function build_impl!( container, transmission, transmission_model.subnetworks, - sys, - ) + sys,) # Order is required for device_model in values(template.devices) From 2cb0c3611842cece3566e3b7dace86195207ebeb Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 28 Dec 2023 17:41:01 -0700 Subject: [PATCH 20/59] add bus maps --- src/core/network_model.jl | 2 ++ src/core/optimization_container.jl | 28 +++++++++++++++++++++++----- 2 files changed, 25 insertions(+), 5 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 32a6f2be32..8c5a21abae 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,6 +33,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches + bus_aggregation_map::Dict{Int, Int} reduce_radial_branches::Bool function NetworkModel( @@ -51,6 +52,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), + Dict{Int, Int}(), reduce_radial_branches, ) end diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 89ba392b50..82e5fb88b2 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -431,9 +431,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractPowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) - ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.ACBus) => _make_container_array(ac_bus_numbers, time_steps), @@ -450,9 +455,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractActivePowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) - ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.ACBus) => _make_container_array(ac_bus_numbers, time_steps), @@ -467,6 +477,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, ::Vector{Int}, ::Type{CopperPlatePowerModel}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) time_steps = get_time_steps(container) subnetworks_ref_buses = collect(keys(subnetworks)) @@ -482,9 +493,14 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{T}, + bus_reduction_map::Dict{Int64, Set{Int64}} ) where {T <: Union{PTDFPowerModel, StandardPTDFModel}} time_steps = get_time_steps(container) - ac_bus_numbers = sort!(collect(Iterators.flatten(values(subnetworks)))) + if isempty(bus_reduction_map) + ac_bus_numbers = collect(Iterators.flatten(values(subnetworks))) + else + ac_bus_numbers = collect(keys(bus_reduction_map)) + end subnetworks = collect(keys(subnetworks)) container.expressions = Dict( ExpressionKey(ActivePowerBalance, PSY.System) => @@ -502,9 +518,10 @@ function initialize_system_expressions!( ::Type{T}, subnetworks::Dict{Int, Set{Int}}, system::PSY.System, + bus_reduction_map::Dict{Int64, Set{Int64}} ) where {T <: PM.AbstractPowerModel} dc_bus_numbers = [PSY.get_number(b) for b in PSY.get_components(PSY.DCBus, system)] - _make_system_expressions!(container, subnetworks, dc_bus_numbers, T) + _make_system_expressions!(container, subnetworks, dc_bus_numbers, T, bus_reduction_map) return end @@ -519,7 +536,8 @@ function build_impl!( container, transmission, transmission_model.subnetworks, - sys,) + sys, + transmission_model.radial_branches.bus_reduction_map) # Order is required for device_model in values(template.devices) From 869042fdb9f15e0a94768509b3bc03c91d10e27d Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 08:34:03 -0700 Subject: [PATCH 21/59] add missing sort! call --- src/core/optimization_container.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 82e5fb88b2..f6074b8e64 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -508,7 +508,7 @@ function _make_system_expressions!( ExpressionKey(ActivePowerBalance, PSY.DCBus) => _make_container_array(dc_bus_numbers, time_steps), ExpressionKey(ActivePowerBalance, PSY.ACBus) => - _make_container_array(ac_bus_numbers, time_steps), + _make_container_array(sort!(ac_bus_numbers), time_steps), ) return end From 72ae721ed420bae8fb2b7c76fb66cbebdfdac4c8 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:10 -0700 Subject: [PATCH 22/59] updates to network model --- src/core/network_model.jl | 5 +++-- src/initial_conditions/initialization.jl | 4 ++++ 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 8c5a21abae..15c0edbee6 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -33,7 +33,6 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} bus_area_map::Dict{PSY.ACBus, Int} duals::Vector{DataType} radial_branches::PNM.RadialBranches - bus_aggregation_map::Dict{Int, Int} reduce_radial_branches::Bool function NetworkModel( @@ -52,7 +51,6 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} Dict{PSY.ACBus, Int}(), duals, PNM.RadialBranches(), - Dict{Int, Int}(), reduce_radial_branches, ) end @@ -60,6 +58,8 @@ 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_duals(m::NetworkModel) = m.duals get_network_formulation(::NetworkModel{T}) where {T} = T get_reference_buses(m::NetworkModel{T}) where {T <: PM.AbstractPowerModel} = @@ -116,6 +116,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: _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 diff --git a/src/initial_conditions/initialization.jl b/src/initial_conditions/initialization.jl index 80ff0f1fd7..d866993118 100644 --- a/src/initial_conditions/initialization.jl +++ b/src/initial_conditions/initialization.jl @@ -5,7 +5,11 @@ function get_initial_conditions_template(model::OperationModel) get_network_formulation(model.template); use_slacks = get_use_slacks(get_network_model(model.template)), PTDF_matrix = get_PTDF_matrix(get_network_model(model.template)), + reduce_radial_branches = get_reduce_radial_branches( + get_network_model(model.template), + ), ) + network_model.radial_branches = get_radial_branches(get_network_model(model.template)) network_model.subnetworks = get_subnetworks(get_network_model(model.template)) bus_area_map = get_bus_area_map(get_network_model(model.template)) if !isempty(bus_area_map) From 4e9d04615fecf522afefbbf42893451c4e319de7 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:25 -0700 Subject: [PATCH 23/59] updates to the container --- src/core/optimization_container.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index f6074b8e64..5e68f6e067 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -431,7 +431,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractPowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -455,7 +455,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{<:PM.AbstractActivePowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -477,7 +477,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, ::Vector{Int}, ::Type{CopperPlatePowerModel}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) time_steps = get_time_steps(container) subnetworks_ref_buses = collect(keys(subnetworks)) @@ -493,7 +493,7 @@ function _make_system_expressions!( subnetworks::Dict{Int, Set{Int}}, dc_bus_numbers::Vector{Int}, ::Type{T}, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) where {T <: Union{PTDFPowerModel, StandardPTDFModel}} time_steps = get_time_steps(container) if isempty(bus_reduction_map) @@ -518,7 +518,7 @@ function initialize_system_expressions!( ::Type{T}, subnetworks::Dict{Int, Set{Int}}, system::PSY.System, - bus_reduction_map::Dict{Int64, Set{Int64}} + bus_reduction_map::Dict{Int64, Set{Int64}}, ) where {T <: PM.AbstractPowerModel} dc_bus_numbers = [PSY.get_number(b) for b in PSY.get_components(PSY.DCBus, system)] _make_system_expressions!(container, subnetworks, dc_bus_numbers, T, bus_reduction_map) From b3d7511ab45ee6a7713ca9c2f3cb8e7b572379f7 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 14:52:33 -0700 Subject: [PATCH 24/59] updates to branch modeling --- src/devices_models/devices/AC_branches.jl | 49 ++++++++- .../devices/common/add_to_expression.jl | 99 +++++++++++++------ 2 files changed, 112 insertions(+), 36 deletions(-) diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 1176f83700..bf24ae3c3c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -168,14 +168,55 @@ function add_constraints!( container::OptimizationContainer, cons_type::Type{RateLimitConstraint}, devices::IS.FlattenIteratorWrapper{T}, - model::DeviceModel{T, U}, - ::NetworkModel{X}, + ::DeviceModel{T, U}, + network_model::NetworkModel{V}, ) where { T <: PSY.ACBranch, U <: AbstractBranchFormulation, - X <: PM.AbstractActivePowerModel, + V <: PM.AbstractActivePowerModel, } - add_range_constraints!(container, cons_type, FlowActivePowerVariable, devices, model, X) + time_steps = get_time_steps(container) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + device_names = [PSY.get_name(d) for d in devices] + else + device_names = PNM.get_meshed_branches(radial_branches) + end + + con_lb = + add_constraints_container!( + container, + cons_type(), + T, + device_names, + time_steps; + meta = "lb", + ) + con_ub = + add_constraints_container!( + container, + cons_type(), + T, + device_names, + time_steps; + meta = "ub", + ) + + array = get_variable(container, FlowActivePowerVariable(), T) + + for device in devices + ci_name = PSY.get_name(device) + if ci_name ∈ PNM.get_radial_branches(radial_branches) + continue + end + limits = get_min_max_limits(device, RateLimitConstraint, U) # depends on constraint type and formulation type + for t in time_steps + con_ub[ci_name, t] = + JuMP.@constraint(container.JuMPmodel, array[ci_name, t] <= limits.max) + con_lb[ci_name, t] = + JuMP.@constraint(container.JuMPmodel, array[ci_name, t] >= limits.min) + end + end return end diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 54aa14a590..22b56059b2 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -79,7 +79,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, model::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: TimeSeriesParameter, @@ -89,11 +89,12 @@ function add_to_expression!( } param_container = get_parameter(container, U(), V) multiplier = get_multiplier_array(param_container) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) name = PSY.get_name(d) _add_to_jump_expression!( - get_expression(container, T(), PSY.ACBus)[bus_number, t], + get_expression(container, T(), PSY.ACBus)[bus_no, t], get_parameter_column_refs(param_container, name)[t], multiplier[name, t], ) @@ -107,7 +108,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: ActivePowerBalance, U <: OnStatusParameter, @@ -116,13 +117,13 @@ function add_to_expression!( X <: PM.AbstractPowerModel, } parameter = get_parameter_array(container, U(), V) - + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) name = PSY.get_name(d) mult = get_expression_multiplier(U(), T(), d, W()) _add_to_jump_expression!( - get_expression(container, T(), PSY.ACBus)[bus_number, t], + get_expression(container, T(), PSY.ACBus)[bus_no, t], parameter[name, t], mult, ) @@ -139,7 +140,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: VariableType, @@ -149,11 +150,12 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, PSY.get_number(PSY.get_bus(d))) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -215,8 +217,10 @@ 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) for d in devices - bus_no_to = PSY.get_number(PSY.get_arc(d).to) + bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) + bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_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) @@ -250,8 +254,10 @@ 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) for d in devices - bus_no_from = PSY.get_number(PSY.get_arc(d).from) + bus_no_from_ = PSY.get_number(PSY.get_arc(d).from) + bus_no_from = PNM.get_mapped_bus_number(radial_branches, bus_no_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) @@ -352,12 +358,14 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_arc(d).from) + bus_no_ = PSY.get_number(PSY.get_arc(d).from) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -385,12 +393,14 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_arc(d).to) + bus_no_ = PSY.get_number(PSY.get_arc(d).to) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) for t in get_time_steps(container) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), V, W()), ) @@ -415,11 +425,13 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(network_model) for d in devices, t in get_time_steps(container) name = PSY.get_name(d) - bus_number = PSY.get_number(PSY.get_bus(d)) + bus_no_ = PSY.get_number(PSY.get_bus(d)) + PNM.get_mapped_bus_number(radial_branches, bus_no_) _add_to_jump_expression!( - expression[bus_number, t], + expression[bus_no, t], variable[name, t], get_variable_multiplier(U(), d, W()), ) @@ -574,10 +586,12 @@ function add_to_expression!( 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) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) param = get_parameter_column_refs(param_container, name) for t in get_time_steps(container) @@ -605,9 +619,11 @@ function add_to_expression!( 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) 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_ = PSY.get_number(PSY.get_bus(d)) + bus_no = PNM.get_mapped_bus_number(radial_branches, 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) @@ -637,10 +653,13 @@ function add_to_expression!( 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) + @assert !isempty(radial_branches) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -675,10 +694,12 @@ function add_to_expression!( 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) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) - bus_no = PSY.get_number(device_bus) + 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -715,16 +736,21 @@ function add_to_expression!( } var = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) + radial_branches = get_radial_branches(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_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).from), t], + expression[bus_no_from, t], flow_variable, -1.0, ) _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).to), t], + expression[bus_no_to, t], flow_variable, 1.0, ) @@ -752,9 +778,12 @@ 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) for d in devices - bus_no_from = PSY.get_number(PSY.get_arc(d).from) - bus_no_to = PSY.get_number(PSY.get_arc(d).to) + 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_) 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) @@ -826,16 +855,21 @@ function add_to_expression!( ) 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) 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_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).from), t], + expression[bus_no_from, t], flow_variable, -get_variable_multiplier(U(), d, V()), ) _add_to_jump_expression!( - expression[PSY.get_number(PSY.get_arc(d).to), t], + expression[bus_no_to, t], flow_variable, get_variable_multiplier(U(), d, V()), ) @@ -1079,6 +1113,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus) expression = get_expression(container, T(), PSY.ACBus) @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(variable, 1) _add_to_jump_expression!( expression[n, t], @@ -1102,8 +1137,8 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus, "P") expression = get_expression(container, T(), PSY.ACBus) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) - for t in get_time_steps(container), n in bus_numbers + # 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(variable, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1126,8 +1161,8 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus, "Q") expression = get_expression(container, T(), PSY.ACBus) - bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) - for t in get_time_steps(container), n in bus_numbers + # 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(variable, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], From 8e1a85c9efd160690623e2ca1d96bd6f8dce8d8b Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 16:32:53 -0700 Subject: [PATCH 25/59] remove debugging assertion --- src/devices_models/devices/common/add_to_expression.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 22b56059b2..f37e4e729c 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -654,7 +654,6 @@ function add_to_expression!( sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) radial_branches = get_radial_branches(network_model) - @assert !isempty(radial_branches) for d in devices name = PSY.get_name(d) device_bus = PSY.get_bus(d) From e1aab9e5bbe65adfe8b24d1018af07785a7df8ee Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 17:11:16 -0700 Subject: [PATCH 26/59] fix tests --- src/devices_models/devices/common/add_to_expression.jl | 6 +++--- test/test_utils/mock_operation_models.jl | 1 + 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index f37e4e729c..31340271dd 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -383,7 +383,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: ActivePowerBalance, U <: FlowActivePowerToFromVariable, @@ -415,7 +415,7 @@ function add_to_expression!( ::Type{U}, devices::IS.FlattenIteratorWrapper{V}, ::DeviceModel{V, W}, - ::NetworkModel{X}, + network_model::NetworkModel{X}, ) where { T <: SystemBalanceExpressions, U <: OnVariable, @@ -429,7 +429,7 @@ function add_to_expression!( for d in devices, t in get_time_steps(container) name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) - PNM.get_mapped_bus_number(radial_branches, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_branches, bus_no_) _add_to_jump_expression!( expression[bus_no, t], variable[name, t], diff --git a/test/test_utils/mock_operation_models.jl b/test/test_utils/mock_operation_models.jl index a41128efb5..759684109c 100644 --- a/test/test_utils/mock_operation_models.jl +++ b/test/test_utils/mock_operation_models.jl @@ -117,6 +117,7 @@ function mock_construct_device!( PSI.get_network_formulation(template), PSI.get_network_model(template).subnetworks, PSI.get_system(problem), + Dict{Int64, Set{Int64}}(), ) if PSI.validate_available_devices(model, PSI.get_system(problem)) PSI.construct_device!( From 55704f10faadf292508dc38c2c10ad2d58daf44f Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 2 Jan 2024 21:26:19 -0700 Subject: [PATCH 27/59] fix bug with serialize function --- src/operation/operation_model_interface.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/operation/operation_model_interface.jl b/src/operation/operation_model_interface.jl index 9a1e01c0c0..faa323ae6b 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -415,7 +415,7 @@ end function serialize_optimization_model(model::OperationModel, save_path::String) serialize_jump_optimization_model( - get_optimization_container(model), + get_jump_model(get_optimization_container(model)), save_path, ) return From 67cf974342de3cb40a6cc62a7c1f4a433a168ca1 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 13:40:20 -0800 Subject: [PATCH 28/59] update branch constructors to ignore radial branches --- .../device_constructors/branch_constructor.jl | 8 ++-- src/devices_models/devices/AC_branches.jl | 37 ++++++++++++++++--- src/network_models/pm_translator.jl | 16 +++++--- src/network_models/powermodels_interface.jl | 35 +++++++++++++----- 4 files changed, 73 insertions(+), 23 deletions(-) diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index ec1b7d6d9d..616500cafb 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -106,7 +106,7 @@ function construct_device!( container, devices, device_model, - get_network_formulation(network_model), + network_model, ) end return @@ -273,7 +273,7 @@ function construct_device!( container, devices, device_model, - get_network_formulation(network_model), + network_model, ) add_constraint_dual!(container, sys, device_model) return @@ -331,7 +331,7 @@ function construct_device!( ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) - branch_rate_bounds!(container, devices, model, get_network_formulation(network_model)) + branch_rate_bounds!(container, devices, model, network_model) add_constraints!(container, RateLimitConstraintFromTo, devices, model, network_model) add_constraints!(container, RateLimitConstraintToFrom, devices, model, network_model) @@ -356,7 +356,7 @@ function construct_device!( ) where {T <: PSY.ACBranch} devices = get_available_components(T, sys, get_attribute(model, "filter_function")) - branch_rate_bounds!(container, devices, model, get_network_formulation(network_model)) + branch_rate_bounds!(container, devices, model, network_model) add_constraint_dual!(container, sys, model) return end diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index bf24ae3c3c..de74ab915c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -104,11 +104,18 @@ function branch_rate_bounds!( container::OptimizationContainer, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::Type{<:PM.AbstractDCPModel}, + network_model::NetworkModel{<:PM.AbstractDCPModel}, ) 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) + for d in devices name = PSY.get_name(d) + if name ∈ radial_branches_names + continue + end for t in get_time_steps(container) JuMP.set_upper_bound(var[name, t], PSY.get_rate(d)) JuMP.set_lower_bound(var[name, t], -1.0 * PSY.get_rate(d)) @@ -121,15 +128,23 @@ function branch_rate_bounds!( container::OptimizationContainer, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::Type{<:PM.AbstractPowerModel}, + network_model::NetworkModel{<:PM.AbstractPowerModel}, ) where {B <: PSY.ACBranch} vars = [ get_variable(container, FlowActivePowerFromToVariable(), B), get_variable(container, FlowActivePowerToFromVariable(), B), ] + + time_steps = get_time_steps(container) + radial_branches = get_radial_branches(network_model) + radial_branches_names = PNM.get_radial_branches(radial_branches) + for d in devices name = PSY.get_name(d) - for t in get_time_steps(container), var in vars + if name ∈ radial_branches_names + continue + end + for t in time_steps, var in vars JuMP.set_upper_bound(var[name, t], PSY.get_rate(d)) JuMP.set_lower_bound(var[name, t], -1.0 * PSY.get_rate(d)) end @@ -256,7 +271,7 @@ function add_constraints!( cons_type::Type{RateLimitConstraintFromTo}, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::NetworkModel{T}, + network_model::NetworkModel{T}, ) where {B <: PSY.ACBranch, T <: PM.AbstractPowerModel} rating_data = [(PSY.get_name(h), PSY.get_rate(h)) for h in devices] @@ -272,7 +287,13 @@ 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) + for r in rating_data + if r[1] ∈ radial_branches_names + continue + end for t in time_steps constraint[r[1], t] = JuMP.@constraint( get_jump_model(container), @@ -291,7 +312,7 @@ function add_constraints!( cons_type::Type{RateLimitConstraintToFrom}, devices::IS.FlattenIteratorWrapper{B}, ::DeviceModel{B, <:AbstractBranchFormulation}, - ::NetworkModel{T}, + network_model::NetworkModel{T}, ) where {B <: PSY.ACBranch, T <: PM.AbstractPowerModel} rating_data = [(PSY.get_name(h), PSY.get_rate(h)) for h in devices] @@ -307,7 +328,13 @@ 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) + for r in rating_data + if r[1] ∈ radial_branches_names + continue + end for t in time_steps constraint[r[1], t] = JuMP.@constraint( get_jump_model(container), diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 627a5b7cf0..ec012ff097 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -381,7 +381,7 @@ end function get_branches_to_pm( sys::PSY.System, - ::Type{S}, + network_model::NetworkModel{S}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -389,6 +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) for (d, device_model) in branch_template comp_type = get_component_type(device_model) if comp_type <: TwoTerminalHVDCTypes @@ -398,6 +400,10 @@ function get_branches_to_pm( start_idx += length(PM_branches) 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" + continue + end ix = i + start_idx PM_branches["$(ix)"] = get_branch_to_pm(ix, branch, get_formulation(device_model), S) @@ -413,7 +419,7 @@ end function get_branches_to_pm( sys::PSY.System, - ::Type{S}, + network_model::NetworkModel{S}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -442,7 +448,7 @@ end function get_branches_to_pm( ::PSY.System, - ::Type{PTDFPowerModel}, + network_model::NetworkModel{PTDFPowerModel}, ::Type{T}, branch_template::BranchModelContainer, start_idx = 0, @@ -485,13 +491,13 @@ end function pass_to_pm(sys::PSY.System, template::ProblemTemplate, time_periods::Int) ac_lines, PMmap_ac = get_branches_to_pm( sys, - get_network_formulation(template), + get_network_model(template), PSY.ACBranch, template.branches, ) two_terminal_dc_lines, PMmap_dc = get_branches_to_pm( sys, - get_network_formulation(template), + get_network_model(template), TwoTerminalHVDCTypes, template.branches, length(ac_lines), diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index 93edc742a5..ff0b9f7ecd 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -290,17 +290,25 @@ function powermodels_network!( ) where {S <: PM.AbstractPowerModel} 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) - for t in time_steps, bus in buses - pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_p"] = + network_model = get_network_model(template) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_reduction_map = radial_branches.bus_reduction_map + ac_bus_numbers = collect(keys(bus_reduction_map)) + end + + for t in time_steps, bus_no in ac_bus_numbers + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_p"] = container.expressions[ExpressionKey(ActivePowerBalance, PSY.ACBus)][ - bus.number, + bus_no, t, ] - pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_q"] = + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_q"] = container.expressions[ExpressionKey(ReactivePowerBalance, PSY.ACBus)][ - bus.number, + bus_no, t, ] end @@ -323,10 +331,19 @@ function powermodels_network!( pm_data, PM_map = pass_to_pm(sys, template, time_steps[end]) buses = get_available_components(PSY.ACBus, sys) - for t in time_steps, bus in buses - pm_data["nw"]["$(t)"]["bus"]["$(PSY.get_number(bus))"]["inj_p"] = + network_model = get_network_model(template) + radial_branches = get_radial_branches(network_model) + if isempty(radial_branches) + ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) + else + bus_reduction_map = radial_branches.bus_reduction_map + ac_bus_numbers = collect(keys(bus_reduction_map)) + end + + for t in time_steps, bus_no in ac_bus_numbers + pm_data["nw"]["$(t)"]["bus"]["$bus_no"]["inj_p"] = container.expressions[ExpressionKey(ActivePowerBalance, PSY.ACBus)][ - PSY.get_number(bus), + bus_no, t, ] # pm_data["nw"]["$(t)"]["bus"]["$(bus.number)"]["inj_q"] = 0.0 From b9195ea1184b18a7e13cc4e2e2b3a247d8265e7a Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 13:40:40 -0800 Subject: [PATCH 29/59] update slack expressions --- src/devices_models/devices/common/add_to_expression.jl | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 31340271dd..17bdefdbb0 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -1111,9 +1111,10 @@ function add_to_expression!( } variable = get_variable(container, U(), PSY.ACBus) expression = get_expression(container, T(), PSY.ACBus) - @assert_op length(axes(variable, 1)) == length(axes(expression, 1)) + #TODO: Update Slacks to consider reduced buses + #@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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1137,7 +1138,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus, "P") expression = get_expression(container, T(), PSY.ACBus) # 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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], @@ -1161,7 +1162,7 @@ function add_to_expression!( variable = get_variable(container, U(), PSY.ACBus, "Q") expression = get_expression(container, T(), PSY.ACBus) # 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(variable, 1) + for t in get_time_steps(container), n in axes(expression, 1) _add_to_jump_expression!( expression[n, t], variable[n, t], From 1901f14d1c651ee442b372ebf449b64cdc8d0a24 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 3 Jan 2024 15:28:08 -0800 Subject: [PATCH 30/59] add radial test --- test/test_radial_branches.jl | 234 +++++++++++++++++++++++++++++++++++ 1 file changed, 234 insertions(+) create mode 100644 test/test_radial_branches.jl diff --git a/test/test_radial_branches.jl b/test/test_radial_branches.jl new file mode 100644 index 0000000000..1ba2633443 --- /dev/null +++ b/test/test_radial_branches.jl @@ -0,0 +1,234 @@ +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 + new_sys = _updated_5bus_sys_with_extensions() + + 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) + + 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 From 46122d53ed33c2714249eff931d8e578cd0b5e1d Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 8 Jan 2024 10:12:22 -0700 Subject: [PATCH 31/59] add new test --- test/test_network_constructors.jl | 62 ++++++++++++++++++++----------- 1 file changed, 40 insertions(+), 22 deletions(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 6e4860a9d3..b7f7fdb8b8 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -1,29 +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), +] + @testset "All PowerModels models construction" begin - 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), - ] c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") - for (network, solver) in networks + for (network, solver) in networks_for_testing template = get_thermal_dispatch_template_network( NetworkModel(network; PTDF_matrix = PTDF(c_sys5)), ) @@ -696,3 +697,20 @@ end ) @test all(isapprox.(sum(zone_2_gen .+ zone_2_load; dims = 2), 0.0; atol = 1e-3)) end + +@testset "All PowerModels models construction with reduced radial branches" begin + # use the system with a radial branch + #c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") + for (network, solver) in networks_for_testing + template = get_thermal_dispatch_template_network( + NetworkModel(network; + PTDF_matrix = PTDF(c_sys5), + reduce_radial_branches = true, + use_slacks), + ) + ps_model = DecisionModel(template, c_sys5; optimizer = solver) + @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == + PSI.BuildStatus.BUILT + @test ps_model.internal.container.pm !== nothing + end +end From f2a2a8859e08142adeadb712e63f230fcaf81bdf Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 8 Jan 2024 10:13:24 -0700 Subject: [PATCH 32/59] update testing --- test/test_network_constructors.jl | 237 +++++++++++++++++++++++++++++- test/test_radial_branches.jl | 234 ----------------------------- 2 files changed, 235 insertions(+), 236 deletions(-) delete mode 100644 test/test_radial_branches.jl diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index b7f7fdb8b8..8a23dc591f 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -698,9 +698,242 @@ 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) + + 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 "All PowerModels models construction with reduced radial branches" begin - # use the system with a radial branch - #c_sys5 = PSB.build_system(PSITestSystems, "c_sys5") + c_sys5 = _updated_5bus_sys_with_extensions() for (network, solver) in networks_for_testing template = get_thermal_dispatch_template_network( NetworkModel(network; diff --git a/test/test_radial_branches.jl b/test/test_radial_branches.jl deleted file mode 100644 index 1ba2633443..0000000000 --- a/test/test_radial_branches.jl +++ /dev/null @@ -1,234 +0,0 @@ -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 - new_sys = _updated_5bus_sys_with_extensions() - - 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) - - 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 From 6ea972435e33ec2b0493cca2ac767ac68dd732b5 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 11 Jan 2024 10:15:55 -0700 Subject: [PATCH 33/59] bump PNM version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 85dbbb0c35..3f7314a724 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,7 @@ LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" PowerModels = "~0.20" -PowerNetworkMatrices = "^0.9" +PowerNetworkMatrices = "^0.10" PowerSystems = "^3" PrettyTables = "2" ProgressMeter = "^1.5" From 59831034079e3a3bce89ab5cbeb12f4a51579ecc Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 11 Jan 2024 16:37:11 -0700 Subject: [PATCH 34/59] bump ts version --- Project.toml | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Project.toml b/Project.toml index 3f7314a724..60d92626ca 100644 --- a/Project.toml +++ b/Project.toml @@ -45,13 +45,13 @@ JuMP = "1" LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" -PowerModels = "~0.20" -PowerNetworkMatrices = "^0.10" +PowerModels = "^0.20" +PowerNetworkMatrices = "^0.9, ^0.10" 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" From 74a643c2f5a8d8c8906f2d03954083bceeefeb33 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 12:31:39 -0700 Subject: [PATCH 35/59] formatter and test fixes --- src/network_models/network_constructor.jl | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/src/network_models/network_constructor.jl b/src/network_models/network_constructor.jl index e2d4a5d850..7f53454439 100644 --- a/src/network_models/network_constructor.jl +++ b/src/network_models/network_constructor.jl @@ -202,15 +202,14 @@ function construct_network!( end if get_use_slacks(model) - add_variables!(container, SystemBalanceSlackUp, sys, T) - add_variables!(container, SystemBalanceSlackDown, sys, T) + add_variables!(container, SystemBalanceSlackUp, sys, model) + add_variables!(container, SystemBalanceSlackDown, sys, model) add_to_expression!( container, ActivePowerBalance, SystemBalanceSlackUp, sys, model, - T, ) add_to_expression!( container, @@ -218,7 +217,6 @@ function construct_network!( SystemBalanceSlackDown, sys, model, - T, ) add_to_expression!( container, @@ -226,7 +224,6 @@ function construct_network!( SystemBalanceSlackUp, sys, model, - T, ) add_to_expression!( container, @@ -234,7 +231,6 @@ function construct_network!( SystemBalanceSlackDown, sys, model, - T, ) objective_function!(container, PSY.ACBus, model) end From d00ad503cac55793c5ff32272b3093a9a5df6cab Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 12:31:51 -0700 Subject: [PATCH 36/59] update te testing code --- test/test_network_constructors.jl | 124 ++---------------------------- 1 file changed, 8 insertions(+), 116 deletions(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 8a23dc591f..6de58ca450 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -1,6 +1,6 @@ # 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 = [ +const networks_for_testing = [ (PM.ACPPowerModel, fast_ipopt_optimizer), (PM.ACRPowerModel, fast_ipopt_optimizer), (PM.ACTPowerModel, fast_ipopt_optimizer), @@ -698,117 +698,8 @@ 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() + new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") net_model = StandardPTDFModel @@ -871,7 +762,7 @@ end end @testset "DCPPowerModel Radial Branches Test" begin - + new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") net_model = DCPPowerModel template_uc = template_unit_commitment(; @@ -933,15 +824,16 @@ end end @testset "All PowerModels models construction with reduced radial branches" begin - c_sys5 = _updated_5bus_sys_with_extensions() + new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") for (network, solver) in networks_for_testing + @error network solver template = get_thermal_dispatch_template_network( NetworkModel(network; - PTDF_matrix = PTDF(c_sys5), + PTDF_matrix = PTDF(new_sys), reduce_radial_branches = true, - use_slacks), + 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 From 32b40cb032d0f263831d7edd7fa41638eb454777 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 12:58:05 -0700 Subject: [PATCH 37/59] rename PNM object --- src/core/network_model.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 15c0edbee6..e11336a7fd 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -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_branches::PNM.RadialNetworkReduction reduce_radial_branches::Bool function NetworkModel( @@ -50,7 +50,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} subnetworks, Dict{PSY.ACBus, Int}(), duals, - PNM.RadialBranches(), + PNM.RadialNetworkReduction(), reduce_radial_branches, ) end @@ -83,7 +83,7 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches - model.radial_branches = PNM.RadialBranches(sys) + model.radial_branches = PNM.RadialNetworkReduction(sys) end return end From f0a7c76cf6c6ffd3e191fbdcba7fa31e12806ff5 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:05:18 -0700 Subject: [PATCH 38/59] replace radial_branches to radial_network_reduction --- src/core/network_model.jl | 8 ++--- src/devices_models/devices/AC_branches.jl | 10 +++---- .../devices/common/add_to_expression.jl | 30 +++++++++---------- src/initial_conditions/initialization.jl | 2 +- src/network_models/pm_translator.jl | 2 +- src/network_models/powermodels_interface.jl | 4 +-- 6 files changed, 28 insertions(+), 28 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index e11336a7fd..1bbb999fbb 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -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.RadialNetworkReduction + radial_network_reduction::PNM.RadialNetworkReduction reduce_radial_branches::Bool function NetworkModel( @@ -59,7 +59,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_branches get_duals(m::NetworkModel) = m.duals get_network_formulation(::NetworkModel{T}) where {T} = T get_reference_buses(m::NetworkModel{T}) where {T <: PM.AbstractPowerModel} = @@ -83,7 +83,7 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches - model.radial_branches = PNM.RadialNetworkReduction(sys) + model.radial_network_reduction = PNM.RadialNetworkReduction(sys) end return end @@ -117,7 +117,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: end if model.reduce_radial_branches @assert !isempty(model.PTDF_matrix.radial_branches) - model.radial_branches = model.PTDF_matrix.radial_branches + model.radial_network_reduction = model.PTDF_matrix.radial_branches end return end diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index de74ab915c..d4e2dd3a17 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -108,7 +108,7 @@ function branch_rate_bounds!( ) where {B <: PSY.ACBranch} var = get_variable(container, FlowActivePowerVariable(), B) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_branches) for d in devices @@ -136,7 +136,7 @@ function branch_rate_bounds!( ] time_steps = get_time_steps(container) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_branches) for d in devices @@ -191,7 +191,7 @@ function add_constraints!( V <: PM.AbstractActivePowerModel, } time_steps = get_time_steps(container) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) if isempty(radial_branches) device_names = [PSY.get_name(d) for d in devices] else @@ -287,7 +287,7 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_branches) for r in rating_data @@ -328,7 +328,7 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_branches) for r in rating_data diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 17bdefdbb0..8a5f5b843d 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -89,7 +89,7 @@ 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))) name = PSY.get_name(d) @@ -117,7 +117,7 @@ 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))) name = PSY.get_name(d) @@ -150,7 +150,7 @@ 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))) @@ -217,7 +217,7 @@ 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_to_ = PSY.get_number(PSY.get_arc(d).to) bus_no_to = PNM.get_mapped_bus_number(radial_branches, bus_no_to_) @@ -254,7 +254,7 @@ 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_) @@ -358,7 +358,7 @@ 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) @@ -393,7 +393,7 @@ 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) @@ -425,7 +425,7 @@ 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)) @@ -586,7 +586,7 @@ function add_to_expression!( 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) @@ -619,7 +619,7 @@ function add_to_expression!( 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)) @@ -653,7 +653,7 @@ function add_to_expression!( 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) @@ -693,7 +693,7 @@ function add_to_expression!( 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) @@ -735,7 +735,7 @@ 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_) @@ -777,7 +777,7 @@ 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_) @@ -854,7 +854,7 @@ function add_to_expression!( ) 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_) diff --git a/src/initial_conditions/initialization.jl b/src/initial_conditions/initialization.jl index d866993118..02a3303e40 100644 --- a/src/initial_conditions/initialization.jl +++ b/src/initial_conditions/initialization.jl @@ -9,7 +9,7 @@ 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)) bus_area_map = get_bus_area_map(get_network_model(model.template)) if !isempty(bus_area_map) diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index ec012ff097..78c686be10 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -389,7 +389,7 @@ 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_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_branches) for (d, device_model) in branch_template comp_type = get_component_type(device_model) diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index ff0b9f7ecd..817211c5be 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -292,7 +292,7 @@ function powermodels_network!( 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) + radial_network_reduction = get_radial_network_reduction(network_model) if isempty(radial_branches) ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) else @@ -332,7 +332,7 @@ function powermodels_network!( buses = get_available_components(PSY.ACBus, sys) network_model = get_network_model(template) - radial_branches = get_radial_branches(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) if isempty(radial_branches) ac_bus_numbers = PSY.get_number.(get_available_components(PSY.ACBus, sys)) else From 2f6b9f0f235c0a1c81176b5fa3f3d4eeb9a98792 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:06:32 -0700 Subject: [PATCH 39/59] replace field names --- src/core/network_model.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 1bbb999fbb..149e447f38 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -59,7 +59,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_network_reduction(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} = From 82865ab3945e81076411e49cc943a79e3f124d13 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:07:07 -0700 Subject: [PATCH 40/59] field name change --- src/core/optimization_container.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 5e68f6e067..4e6f2b3fec 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -537,7 +537,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) From a5764c39d6a971157827f865c5e3b6e7700e3acf Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:10:49 -0700 Subject: [PATCH 41/59] more name replace --- src/devices_models/devices/AC_branches.jl | 14 +++---- .../devices/common/add_to_expression.jl | 38 +++++++++---------- src/network_models/pm_translator.jl | 2 +- src/network_models/powermodels_interface.jl | 4 +- 4 files changed, 28 insertions(+), 30 deletions(-) diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index d4e2dd3a17..26821565e6 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -109,7 +109,7 @@ function branch_rate_bounds!( var = get_variable(container, FlowActivePowerVariable(), B) radial_network_reduction = get_radial_network_reduction(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices name = PSY.get_name(d) @@ -137,7 +137,7 @@ function branch_rate_bounds!( time_steps = get_time_steps(container) radial_network_reduction = get_radial_network_reduction(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices name = PSY.get_name(d) @@ -192,10 +192,10 @@ function add_constraints!( } time_steps = get_time_steps(container) radial_network_reduction = get_radial_network_reduction(network_model) - if isempty(radial_branches) + 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 @@ -288,7 +288,7 @@ function add_constraints!( constraint = get_constraint(container, cons_type(), B) radial_network_reduction = get_radial_network_reduction(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data if r[1] ∈ radial_branches_names @@ -329,7 +329,7 @@ function add_constraints!( constraint = get_constraint(container, cons_type(), B) radial_network_reduction = get_radial_network_reduction(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data if r[1] ∈ radial_branches_names diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 8a5f5b843d..efa902cf12 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -91,7 +91,7 @@ function add_to_expression!( multiplier = get_multiplier_array(param_container) 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], @@ -119,7 +119,7 @@ function add_to_expression!( parameter = get_parameter_array(container, U(), V) 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!( @@ -153,7 +153,7 @@ function add_to_expression!( 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], @@ -219,8 +219,7 @@ function add_to_expression!( sys_expr = get_expression(container, T(), PSY.System) 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) @@ -256,8 +255,7 @@ function add_to_expression!( sys_expr = get_expression(container, T(), PSY.System) 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) @@ -362,7 +360,7 @@ function add_to_expression!( 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], @@ -397,7 +395,7 @@ function add_to_expression!( 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], @@ -429,7 +427,7 @@ function add_to_expression!( 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], @@ -591,7 +589,7 @@ function add_to_expression!( 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) @@ -623,7 +621,7 @@ function add_to_expression!( 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) @@ -658,7 +656,7 @@ function add_to_expression!( 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -698,7 +696,7 @@ function add_to_expression!( 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) for t in get_time_steps(container) _add_to_jump_expression!( @@ -738,9 +736,9 @@ function add_to_expression!( 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, 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_to = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_to_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( @@ -780,9 +778,9 @@ function add_to_expression!( 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, 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_to = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_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) @@ -857,9 +855,9 @@ function add_to_expression!( 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, 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_to = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_to_) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] _add_to_jump_expression!( diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 78c686be10..317d5f063a 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -390,7 +390,7 @@ function get_branches_to_pm( PMmap_br = Dict{PM_MAP_TUPLE, T}() radial_network_reduction = get_radial_network_reduction(network_model) - radial_branches_names = PNM.get_radial_branches(radial_branches) + 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 diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index 817211c5be..d15889163b 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -293,7 +293,7 @@ function powermodels_network!( network_model = get_network_model(template) radial_network_reduction = get_radial_network_reduction(network_model) - if isempty(radial_branches) + 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 @@ -333,7 +333,7 @@ function powermodels_network!( network_model = get_network_model(template) radial_network_reduction = get_radial_network_reduction(network_model) - if isempty(radial_branches) + 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 From 07704426b782815ed4b0b185d5a9360ac3f45017 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:16:03 -0700 Subject: [PATCH 42/59] implement name changes --- src/core/network_model.jl | 4 +- src/devices_models/devices/AC_branches.jl | 10 ++-- .../devices/common/add_to_expression.jl | 58 +++++++++---------- src/initial_conditions/initialization.jl | 3 +- src/network_models/pm_translator.jl | 2 +- src/network_models/powermodels_interface.jl | 6 +- test/test_network_constructors.jl | 6 +- 7 files changed, 42 insertions(+), 47 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 149e447f38..bd6747e0d7 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -83,7 +83,7 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches - model.radial_network_reduction = PNM.RadialNetworkReduction(sys) + model.radial_network_reduction = PNM.RadialNetworkReduction(sys) end return end @@ -117,7 +117,7 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: end if model.reduce_radial_branches @assert !isempty(model.PTDF_matrix.radial_branches) - model.radial_network_reduction = model.PTDF_matrix.radial_branches + model.radial_network_reduction = model.PTDF_matrix.radial_branches end return end diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 26821565e6..e14215cec1 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -108,7 +108,7 @@ function branch_rate_bounds!( ) where {B <: PSY.ACBranch} var = get_variable(container, FlowActivePowerVariable(), B) - radial_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices @@ -136,7 +136,7 @@ function branch_rate_bounds!( ] time_steps = get_time_steps(container) - radial_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for d in devices @@ -191,7 +191,7 @@ function add_constraints!( V <: PM.AbstractActivePowerModel, } time_steps = get_time_steps(container) - radial_network_reduction = get_radial_network_reduction(network_model) + 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 @@ -287,7 +287,7 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data @@ -328,7 +328,7 @@ function add_constraints!( ) constraint = get_constraint(container, cons_type(), B) - radial_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) radial_branches_names = PNM.get_radial_branches(radial_network_reduction) for r in rating_data diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index efa902cf12..1fc76185c6 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -89,7 +89,7 @@ function add_to_expression!( } param_container = get_parameter(container, U(), V) multiplier = get_multiplier_array(param_container) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, PSY.get_bus(d)) name = PSY.get_name(d) @@ -117,7 +117,7 @@ function add_to_expression!( X <: PM.AbstractPowerModel, } parameter = get_parameter_array(container, U(), V) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, PSY.get_bus(d)) name = PSY.get_name(d) @@ -150,7 +150,7 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, PSY.get_bus(d)) @@ -217,7 +217,7 @@ 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_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices 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) @@ -253,9 +253,10 @@ 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_network_reduction = get_radial_network_reduction(network_model) + radial_network_reduction = get_radial_network_reduction(network_model) for d in devices - bus_no_from = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_arc(d).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) @@ -356,7 +357,7 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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) @@ -391,7 +392,7 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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) @@ -423,7 +424,7 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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)) @@ -584,7 +585,7 @@ function add_to_expression!( multiplier = get_multiplier_array(param_container) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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) @@ -617,7 +618,7 @@ function add_to_expression!( parameter = get_parameter_array(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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)) @@ -651,7 +652,7 @@ function add_to_expression!( variable = get_variable(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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) @@ -691,12 +692,10 @@ function add_to_expression!( variable = get_variable(container, U(), V) sys_expr = get_expression(container, T(), PSY.System) nodal_expr = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, bus_no_) + bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) ref_bus = get_reference_bus(network_model, device_bus) for t in get_time_steps(container) _add_to_jump_expression!( @@ -733,12 +732,11 @@ function add_to_expression!( } var = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, 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!( @@ -775,12 +773,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_network_reduction = get_radial_network_reduction(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_network_reduction, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, 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) @@ -852,12 +849,11 @@ function add_to_expression!( ) where {T <: ActivePowerBalance, U <: PhaseShifterAngle, V <: PhaseAngleControl} var = get_variable(container, U(), PSY.PhaseShiftingTransformer) expression = get_expression(container, T(), PSY.ACBus) - radial_network_reduction = get_radial_network_reduction(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_network_reduction, bus_no_from_) - bus_no_to_ = PSY.get_number(PSY.get_arc(d).to) - bus_no_to = PNM.get_mapped_bus_number(radial_network_reduction, 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!( diff --git a/src/initial_conditions/initialization.jl b/src/initial_conditions/initialization.jl index 02a3303e40..e5222f6bfc 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_network_reduction = get_radial_network_reduction(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)) bus_area_map = get_bus_area_map(get_network_model(model.template)) if !isempty(bus_area_map) diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 317d5f063a..7037854687 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -389,7 +389,7 @@ function get_branches_to_pm( PM_branches = Dict{String, Any}() PMmap_br = Dict{PM_MAP_TUPLE, T}() - radial_network_reduction = get_radial_network_reduction(network_model) + 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) diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index d15889163b..498641bf21 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -13,7 +13,6 @@ end # replicates PM.build_mn_opf function instantiate_nip_expr(pm::PM.AbstractPowerModel) for n in eachindex(PM.nws(pm)) - 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) @@ -103,7 +102,6 @@ end # replicates PM.build_mn_opf_bf_strg function instantiate_bfp_expr(pm::PM.AbstractPowerModel) for n in eachindex(PM.nws(pm)) - 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) @@ -292,7 +290,7 @@ function powermodels_network!( pm_data, PM_map = pass_to_pm(sys, template, time_steps[end]) network_model = get_network_model(template) - radial_network_reduction = get_radial_network_reduction(network_model) + 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 @@ -332,7 +330,7 @@ function powermodels_network!( buses = get_available_components(PSY.ACBus, sys) network_model = get_network_model(template) - radial_network_reduction = get_radial_network_reduction(network_model) + 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 diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 6de58ca450..e17f47c0ea 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -829,9 +829,9 @@ end @error network solver template = get_thermal_dispatch_template_network( NetworkModel(network; - PTDF_matrix = PTDF(new_sys), - reduce_radial_branches = true, - use_slacks = true), + PTDF_matrix = PTDF(new_sys), + reduce_radial_branches = true, + use_slacks = true), ) ps_model = DecisionModel(template, new_sys; optimizer = solver) @test build!(ps_model; output_dir = mktempdir(; cleanup = true)) == From f856011419bf9f3b7fc988f55052f1424cab8a71 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 13:47:12 -0700 Subject: [PATCH 43/59] fix usage of dot access --- src/network_models/powermodels_interface.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index 498641bf21..f9b6677dc4 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -294,7 +294,7 @@ function powermodels_network!( 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 @@ -334,7 +334,7 @@ function powermodels_network!( 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 From 5a5e9b84766afa86ff0729d1673914fc743340ee Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:22:50 -0700 Subject: [PATCH 44/59] fix incompatible models --- src/core/definitions.jl | 3 --- src/core/network_model.jl | 8 ++++++++ src/devices_models/devices/common/add_to_expression.jl | 3 +-- src/network_models/powermodels_interface.jl | 8 ++++++++ test/test_network_constructors.jl | 4 +++- 5 files changed, 20 insertions(+), 6 deletions(-) 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/network_model.jl b/src/core/network_model.jl index bd6747e0d7..dc22e5e1d7 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -75,6 +75,13 @@ function add_dual!(model::NetworkModel, dual) return end +function check_radial_branch_redution_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, @@ -83,6 +90,7 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches + check_radial_branch_redution_compatibility(T) model.radial_network_reduction = PNM.RadialNetworkReduction(sys) end return diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 1fc76185c6..2c6a3f2726 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -1105,8 +1105,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/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index f9b6677dc4..b999b8b54d 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -6,6 +6,14 @@ ################################################################################# # Model Definitions +const UNSUPPORTED_POWERMODELS = + [PM.SOCBFPowerModel, PM.SOCBFConicPowerModel, PM.IVRPowerModel] + +const INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS = [ + PM.SparseSDPWRMPowerModel, + PTDFPowerModel, + ] + function instantiate_nip_expr_model(data::Dict{String, Any}, model_constructor; kwargs...) return PM.instantiate_model(data, model_constructor, instantiate_nip_expr; kwargs...) end diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index e17f47c0ea..d06b978b5c 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -826,7 +826,9 @@ end @testset "All PowerModels models construction with reduced radial branches" begin new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") for (network, solver) in networks_for_testing - @error network solver + if network ∈ PSI.INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS + continue + end template = get_thermal_dispatch_template_network( NetworkModel(network; PTDF_matrix = PTDF(new_sys), From 4a5a05e53c911ac316872049a132b4d5e2623208 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:26:23 -0700 Subject: [PATCH 45/59] fix typo --- src/network_models/pm_translator.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 7037854687..ee86e0544b 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -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 From 035ef787b24976d0e868847c0e0d207d8d8b391e Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:30:50 -0700 Subject: [PATCH 46/59] use radial network reduction in slacks specs --- src/network_models/network_slack_variables.jl | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/src/network_models/network_slack_variables.jl b/src/network_models/network_slack_variables.jl index 619aeab047..82e1638e49 100644 --- a/src/network_models/network_slack_variables.jl +++ b/src/network_models/network_slack_variables.jl @@ -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 = From 2ab1dd465304777ccb563ef80adc218777478178 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:35:50 -0700 Subject: [PATCH 47/59] remove code duplication --- test/test_network_constructors.jl | 162 ++++++++++-------------------- 1 file changed, 51 insertions(+), 111 deletions(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index d06b978b5c..09af6c517a 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -698,128 +698,68 @@ end @test all(isapprox.(sum(zone_2_gen .+ zone_2_load; dims = 2), 0.0; atol = 1e-3)) end -@testset "StandardPTDF Radial Branches Test" begin +# 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") - net_model = StandardPTDFModel + for net_model in [DCPPowerModel, 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 - new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") - 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) + 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 From c9a7b5a1b64a74508900f81767f3819eeaebdebb Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:38:42 -0700 Subject: [PATCH 48/59] formatter run --- src/core/network_model.jl | 4 +++- src/network_models/powermodels_interface.jl | 6 +++--- test/test_network_constructors.jl | 5 ++--- 3 files changed, 8 insertions(+), 7 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index dc22e5e1d7..d312b56ce9 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -75,7 +75,9 @@ function add_dual!(model::NetworkModel, dual) return end -function check_radial_branch_redution_compatibility(::Type{T}) where {T <: PM.AbstractPowerModel} +function check_radial_branch_redution_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 diff --git a/src/network_models/powermodels_interface.jl b/src/network_models/powermodels_interface.jl index b999b8b54d..ea0b6bb991 100644 --- a/src/network_models/powermodels_interface.jl +++ b/src/network_models/powermodels_interface.jl @@ -10,9 +10,9 @@ const UNSUPPORTED_POWERMODELS = [PM.SOCBFPowerModel, PM.SOCBFConicPowerModel, PM.IVRPowerModel] const INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS = [ - PM.SparseSDPWRMPowerModel, - PTDFPowerModel, - ] + PM.SparseSDPWRMPowerModel, + PTDFPowerModel, +] function instantiate_nip_expr_model(data::Dict{String, Any}, model_constructor; kwargs...) return PM.instantiate_model(data, model_constructor, instantiate_nip_expr; kwargs...) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 09af6c517a..dfa5b4f609 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -703,7 +703,6 @@ end new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") for net_model in [DCPPowerModel, StandardPTDFModel] - template_uc = template_unit_commitment(; network = NetworkModel(net_model; reduce_radial_branches = true, @@ -724,7 +723,7 @@ end ) @test build!(uc_model_red; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT + PSI.BuildStatus.BUILT solve!(uc_model_red) res_red = ProblemResults(uc_model_red) @@ -750,7 +749,7 @@ end ) @test build!(uc_model_orig; output_dir = mktempdir(; cleanup = true)) == - PSI.BuildStatus.BUILT + PSI.BuildStatus.BUILT solve!(uc_model_orig) res_orig = ProblemResults(uc_model_orig) From 3bb6e780e209563d9d6456e10b3c64d3b6c36d84 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 14:40:25 -0700 Subject: [PATCH 49/59] fix eliminated variable --- src/devices_models/devices/common/add_to_expression.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 2c6a3f2726..c3f9e8ba2d 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -696,7 +696,7 @@ function add_to_expression!( for d in devices name = PSY.get_name(d) bus_no = PNM.get_mapped_bus_number(radial_network_reduction, PSY.get_bus(d)) - ref_bus = get_reference_bus(network_model, device_bus) + 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], From 4c1ccbd5b1f509fa0623f8bc7403b302542313c1 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 15 Jan 2024 15:12:35 -0700 Subject: [PATCH 50/59] variable rename --- test/test_network_constructors.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index dfa5b4f609..519b19d28f 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -1,6 +1,6 @@ # 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 = [ +const NETWORKS_FOR_TESTING = [ (PM.ACPPowerModel, fast_ipopt_optimizer), (PM.ACRPowerModel, fast_ipopt_optimizer), (PM.ACTPowerModel, fast_ipopt_optimizer), @@ -24,7 +24,7 @@ const networks_for_testing = [ @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)), ) @@ -764,7 +764,7 @@ end @testset "All PowerModels models construction with reduced radial branches" begin new_sys = PSB.build_system(PSITestSystems, "c_sys5_radial") - for (network, solver) in networks_for_testing + for (network, solver) in NETWORKS_FOR_TESTING if network ∈ PSI.INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS continue end From c92e073b19eccb58b63e49f9dcd2d3f64d1fccdf Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 16 Jan 2024 13:02:44 -0700 Subject: [PATCH 51/59] add missing docstring entry --- src/core/network_model.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index d312b56ce9..290af059aa 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -20,7 +20,7 @@ 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) @@ -35,7 +35,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} radial_network_reduction::PNM.RadialNetworkReduction reduce_radial_branches::Bool - function NetworkModel( + ::Type{T}; use_slacks = false, PTDF_matrix = nothing, From c046be777c5f25ff042b9bb492582cf99383215f Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 16 Jan 2024 13:56:51 -0700 Subject: [PATCH 52/59] fix deleted line --- src/core/network_model.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 290af059aa..ae8634aa66 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -35,7 +35,7 @@ mutable struct NetworkModel{T <: PM.AbstractPowerModel} radial_network_reduction::PNM.RadialNetworkReduction reduce_radial_branches::Bool - + function NetworkModel( ::Type{T}; use_slacks = false, PTDF_matrix = nothing, From db9e04c09e3a46eb1ba2528d06bd95cad28c95b0 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Thu, 18 Jan 2024 13:36:18 -0700 Subject: [PATCH 53/59] rename field --- src/core/network_model.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index ae8634aa66..84b9bdfb2a 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -126,8 +126,8 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: _assign_subnetworks_to_buses(model, sys) end if model.reduce_radial_branches - @assert !isempty(model.PTDF_matrix.radial_branches) - model.radial_network_reduction = model.PTDF_matrix.radial_branches + @assert !isempty(model.PTDF_matrix.radial_network_reduction) + model.radial_network_reduction = model.PTDF_matrix.radial_network_reduction end return end From c6539e3acfb114450d0647589a829acf13fb9e5f Mon Sep 17 00:00:00 2001 From: jd-lara Date: Thu, 18 Jan 2024 20:09:59 -0700 Subject: [PATCH 54/59] use bus --- src/devices_models/devices/common/add_to_expression.jl | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index c3f9e8ba2d..c9c7e1cb68 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -656,8 +656,7 @@ function add_to_expression!( 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_network_reduction, 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!( From d2ba9484cb2ba62d96da79906224c4f85e8904bd Mon Sep 17 00:00:00 2001 From: jd-lara Date: Sun, 21 Jan 2024 21:38:28 -0700 Subject: [PATCH 55/59] fix network instantiation --- src/core/network_model.jl | 20 ++++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index 84b9bdfb2a..cd7fe33bf0 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -119,16 +119,16 @@ function instantiate_network_model(model::NetworkModel{StandardPTDFModel}, sys:: 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_network_reduction) - model.radial_network_reduction = model.PTDF_matrix.radial_network_reduction - end return end @@ -138,19 +138,27 @@ function _assign_subnetworks_to_buses( ) where {T <: Union{CopperPlatePowerModel, StandardPTDFModel}} 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 From e356d1e35980e2e1785d4bb34c3f547b7f99c2a6 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 22 Jan 2024 11:57:36 -0700 Subject: [PATCH 56/59] limit PNM version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 60d92626ca..fd1cf39fb9 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,7 @@ LinearAlgebra = "1" Logging = "1" MathOptInterface = "1" PowerModels = "^0.20" -PowerNetworkMatrices = "^0.9, ^0.10" +PowerNetworkMatrices = "^0.10" PowerSystems = "^3" PrettyTables = "2" ProgressMeter = "^1.5" From 9f1c16d57509bb8d063d360a2c4437b94ea69464 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 22 Jan 2024 13:06:59 -0700 Subject: [PATCH 57/59] fix typo --- src/core/network_model.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/core/network_model.jl b/src/core/network_model.jl index cd7fe33bf0..caa2d8873f 100644 --- a/src/core/network_model.jl +++ b/src/core/network_model.jl @@ -75,7 +75,7 @@ function add_dual!(model::NetworkModel, dual) return end -function check_radial_branch_redution_compatibility( +function check_radial_branch_reduction_compatibility( ::Type{T}, ) where {T <: PM.AbstractPowerModel} if T ∈ INCOMPATIBLE_WITH_RADIAL_BRANCHES_POWERMODELS @@ -92,7 +92,7 @@ function instantiate_network_model( model.subnetworks = PNM.find_subnetworks(sys) end if model.reduce_radial_branches - check_radial_branch_redution_compatibility(T) + check_radial_branch_reduction_compatibility(T) model.radial_network_reduction = PNM.RadialNetworkReduction(sys) end return From 782e0176af5fefc2af18ae9f7f710676283f357c Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 22 Jan 2024 14:02:17 -0700 Subject: [PATCH 58/59] add comment on sorting --- src/core/optimization_container.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 4e6f2b3fec..4611822529 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -508,6 +508,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 From f93f740f8d4ffc61895bf3a8f59f400fdc3961a0 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 22 Jan 2024 14:35:29 -0700 Subject: [PATCH 59/59] Update src/core/optimization_container.jl Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> --- src/core/optimization_container.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index 4611822529..9686a1d06d 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -508,8 +508,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 + # Bus numbers are sorted to guarantee consistency in the order between the + # containers _make_container_array(sort!(ac_bus_numbers), time_steps), ) return