From eb5c88c7e14a18f61a039d3878f7173dc2aa2891 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 7 Aug 2024 17:24:22 -0700 Subject: [PATCH 01/68] add two terminal constraints and variables --- src/core/constraints.jl | 1 + src/core/formulations.jl | 4 + src/core/variables.jl | 7 + .../device_constructors/branch_constructor.jl | 120 ++++++++++++++++++ src/devices_models/devices/AC_branches.jl | 79 ++++++++++++ .../devices/TwoTerminalDC_branches.jl | 56 ++++++-- 6 files changed, 257 insertions(+), 10 deletions(-) diff --git a/src/core/constraints.jl b/src/core/constraints.jl index a36e8257cd..a8300d356d 100644 --- a/src/core/constraints.jl +++ b/src/core/constraints.jl @@ -349,6 +349,7 @@ The specified constraints are formulated as: struct HVDCLossesAbsoluteValue <: ConstraintType end struct HVDCDirection <: ConstraintType end struct InterfaceFlowLimit <: ConstraintType end +struct HVDCFlowCalculationConstraint <: ConstraintType end abstract type PowerVariableLimitsConstraint <: ConstraintType end """ diff --git a/src/core/formulations.jl b/src/core/formulations.jl index 1b6ddb109a..5408d934be 100644 --- a/src/core/formulations.jl +++ b/src/core/formulations.jl @@ -134,6 +134,10 @@ struct HVDCTwoTerminalLossless <: AbstractTwoTerminalDCLineFormulation end Branch type to represent lossy power flow on DC lines """ struct HVDCTwoTerminalDispatch <: AbstractTwoTerminalDCLineFormulation end +""" +Branch type to represent piecewise lossy power flow on two terminal DC lines +""" +struct HVDCTwoTerminalPiecewiseLoss <: AbstractTwoTerminalDCLineFormulation end # Not Implemented # struct VoltageSourceDC <: AbstractTwoTerminalDCLineFormulation end diff --git a/src/core/variables.jl b/src/core/variables.jl index 42d5bc0560..53fe2de1b8 100644 --- a/src/core/variables.jl +++ b/src/core/variables.jl @@ -228,6 +228,13 @@ Docs abbreviation: ``u^\\text{dir}`` """ struct HVDCFlowDirectionVariable <: VariableType end +""" +Struct to dispatch the creation of HVDC Piecewise Loss Variables + +Docs abbreviation: ``h`` +""" +struct HVDCPiecewiseLossVariable <: VariableType end + abstract type SparseVariableType <: VariableType end """ diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index baa4768200..04767cb3f4 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -757,6 +757,126 @@ function construct_device!( return end +function construct_device!( + container::OptimizationContainer, + sys::PSY.System, + ::ArgumentConstructStage, + model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + network_model::NetworkModel{<:AbstractPTDFModel}, +) where {T <: TwoTerminalHVDCTypes} + devices = + get_available_components(model, sys) + add_variables!( + container, + FlowActivePowerToFromVariable, + devices, + HVDCTwoTerminalPiecewiseLoss(), + ) + add_variables!( + container, + FlowActivePowerFromToVariable, + devices, + HVDCTwoTerminalPiecewiseLoss(), + ) + _add_dense_pwl_loss_variables!(container, devices, model) + add_to_expression!( + container, + ActivePowerBalance, + FlowActivePowerToFromVariable, + devices, + model, + network_model, + ) + add_to_expression!( + container, + ActivePowerBalance, + FlowActivePowerFromToVariable, + devices, + model, + network_model, + ) + return +end + +function construct_device!( + container::OptimizationContainer, + sys::PSY.System, + ::ModelConstructStage, + model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + network_model::NetworkModel{<:AbstractPTDFModel}, +) where {T <: TwoTerminalHVDCTypes} + devices = + get_available_components(model, sys) + add_constraints!(container, FlowRateConstraintFromTo, devices, model, network_model) + add_constraints!(container, FlowRateConstraintToFrom, devices, model, network_model) + add_constraints!( + container, + HVDCFlowCalculationConstraint, + devices, + model, + network_model, + ) + return +end + +function construct_device!( + container::OptimizationContainer, + sys::PSY.System, + ::ArgumentConstructStage, + model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + network_model::NetworkModel{<:PM.AbstractActivePowerModel}, +) where {T <: TwoTerminalHVDCTypes} + devices = + get_available_components(model, sys) + add_variables!( + container, + FlowActivePowerToFromVariable, + devices, + HVDCTwoTerminalDispatch(), + ) + add_variables!( + container, + FlowActivePowerFromToVariable, + devices, + HVDCTwoTerminalDispatch(), + ) + add_variables!(container, HVDCFlowDirectionVariable, devices, HVDCTwoTerminalDispatch()) + add_to_expression!( + container, + ActivePowerBalance, + FlowActivePowerToFromVariable, + devices, + model, + network_model, + ) + add_to_expression!( + container, + ActivePowerBalance, + FlowActivePowerFromToVariable, + devices, + model, + network_model, + ) + add_feedforward_arguments!(container, model, devices) + return +end + +function construct_device!( + container::OptimizationContainer, + sys::PSY.System, + ::ModelConstructStage, + model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + network_model::NetworkModel{CopperPlatePowerModel}, +) where {T <: TwoTerminalHVDCTypes} + devices = + get_available_components(model, sys) + @warn "CopperPlatePowerModel models with HVDC ignores inter-area losses" + add_constraints!(container, FlowRateConstraintFromTo, devices, model, network_model) + add_constraints!(container, FlowRateConstraintToFrom, devices, model, network_model) + add_constraint_dual!(container, sys, model) + return +end + function construct_device!( container::OptimizationContainer, sys::PSY.System, diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index fc7b1e1e51..06eb06536c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -32,6 +32,7 @@ get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.ACBranch}, :: # These two methods are defined to avoid ambiguities get_variable_binary(::FlowActivePowerSlackUpperBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false +get_variable_binary(::HVDCPiecewiseLossVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_upper_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = nothing get_variable_lower_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = 0.0 get_variable_upper_bound(::FlowActivePowerSlackLowerBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = nothing @@ -169,6 +170,84 @@ function branch_rate_bounds!( return end +################################## PWL Loss Variables ################################## + +function _check_pwl_loss_model(devices) + first_loss = PSY.get_loss(first(devices)) + first_loss_type = typeof(first_loss) + for d in devices + loss = PSY.get_loss(d) + if !isa(loss, first_loss_type) + error( + "Not all TwoTerminal HVDC lines have the same loss model data. Check that all loss models are LinearCurve or PiecewiseIncrementalCurve", + ) + end + if isa(first_loss, PSY.PiecewiseIncrementalCurve) + len_first_loss = length(PSY.get_slopes(first_loss)) + len_loss = length(PSY.get_slopes(loss)) + if len_first_loss != len_loss + error( + "Different length of PWL segments for TwoTerminal HVDC losses are not supported. Check that all HVDC data have the same amount of PWL segments.", + ) + end + end + end +end + +function _add_dense_pwl_loss_variables!( + container::OptimizationContainer, + devices, + model::DeviceModel{D, HVDCTwoTerminalPiecewiseLoss}, +) where {D <: TwoTerminalHVDCTypes} + # Check if type and length of PWL loss model are the same for all devices + _check_pwl_loss_model(devices) + + # Create Variables + time_steps = get_time_steps(container) + settings = get_settings(container) + formulation = HVDCTwoTerminalPiecewiseLoss() + T = HVDCPiecewiseLossVariable + binary = get_variable_binary(T(), D, formulation) + first_loss = PSY.get_loss(first(devices)) + if isa(first_loss, PSY.LinearCurve) + len_segments = 4 # 2*1 + 2 + elseif isa(first_loss, PSY.PiecewiseIncrementalCurve) + len_segments = 2 * length(PSY.get_slopes(first_loss)) + 2 + else + error("Should not be here") + end + + segments = ["pwl_$i" for i in 1:len_segments] + T = HVDCPiecewiseLossVariable + variable = add_variable_container!( + container, + T(), + D, + [PSY.get_name(d) for d in devices], + segments, + time_steps, + ) + + for t in time_steps, s in segments, d in devices + name = PSY.get_name(d) + variable[name, s, t] = JuMP.@variable( + get_jump_model(container), + base_name = "$(T)_$(D)_{$(name), $(s), $(t)}", + binary = binary + ) + ub = get_variable_upper_bound(T(), d, formulation) + ub !== nothing && JuMP.set_upper_bound(variable[name, s, t], ub) + + lb = get_variable_lower_bound(T(), d, formulation) + lb !== nothing && JuMP.set_lower_bound(variable[name, s, t], lb) + + if get_warm_start(settings) + init = get_variable_warm_start_value(T(), d, formulation) + init !== nothing && JuMP.set_start_value(variable[name, s, t], init) + end + end +end + ################################## Rate Limits constraint_infos ############################ """ diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index b5805456d0..97c693635d 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -42,8 +42,14 @@ function get_variable_multiplier( d::PSY.TwoTerminalHVDCLine, ::HVDCTwoTerminalDispatch, ) - l1 = PSY.get_loss(d).l1 - l0 = PSY.get_loss(d).l0 + loss = PSY.get_loss(d) + if !isa(loss, PSY.LinearCurve) + error( + "HVDCTwoTerminalDispatch of branch $(PSY.get_name(d)) only accepts LinearCurve for loss models.", + ) + end + l1 = PSY.get_proportional_term(loss) + l0 = PSY.get_constant_term(loss) if l1 == 0.0 && l0 == 0.0 return 0.0 else @@ -110,8 +116,14 @@ function get_variable_upper_bound( d::PSY.TwoTerminalHVDCLine, ::HVDCTwoTerminalDispatch, ) - l1 = PSY.get_loss(d).l1 - l0 = PSY.get_loss(d).l0 + loss = PSY.get_loss(d) + if !isa(loss, PSY.LinearCurve) + error( + "HVDCTwoTerminalDispatch of branch $(PSY.get_name(d)) only accepts LinearCurve for loss models.", + ) + end + l1 = PSY.get_proportional_term(loss) + l0 = PSY.get_constant_term(loss) if l1 == 0.0 && l0 == 0.0 return 0.0 else @@ -279,11 +291,12 @@ function add_constraints!( container::OptimizationContainer, ::Type{T}, devices::IS.FlattenIteratorWrapper{U}, - model::DeviceModel{U, HVDCTwoTerminalDispatch}, + model::DeviceModel{U, V}, network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, U <: PSY.TwoTerminalHVDCLine, + V <: Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, } inter_network_branches = U[] for d in devices @@ -303,10 +316,27 @@ function add_constraints!( container::OptimizationContainer, ::Type{T}, devices::IS.FlattenIteratorWrapper{U}, - ::DeviceModel{U, HVDCTwoTerminalDispatch}, + ::DeviceModel{U, V}, ::NetworkModel{<:PM.AbstractDCPModel}, -) where {T <: Union{FlowRateConstraintToFrom, FlowRateConstraintFromTo}, - U <: PSY.TwoTerminalHVDCLine} +) where { + T <: Union{FlowRateConstraintToFrom, FlowRateConstraintFromTo}, + U <: PSY.TwoTerminalHVDCLine, + V <: Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, +} + _add_hvdc_flow_constraints!(container, devices, T()) + return +end + +function add_constraints!( + container::OptimizationContainer, + ::Type{T}, + devices::IS.FlattenIteratorWrapper{U}, + ::DeviceModel{U, HVDCTwoTerminalDispatch}, + ::NetworkModel{<:AbstractPTDFModel}, +) where { + T <: Union{FlowRateConstraintToFrom, FlowRateConstraintFromTo}, + U <: PSY.TwoTerminalHVDCLine, +} _add_hvdc_flow_constraints!(container, devices, T()) return end @@ -405,8 +435,14 @@ function add_constraints!( meta = "ft_lb", ) for d in devices - l1 = PSY.get_loss(d).l1 - l0 = PSY.get_loss(d).l0 + loss = PSY.get_loss(d) + if !isa(loss, PSY.LinearCurve) + error( + "HVDCTwoTerminalDispatch of branch $(PSY.get_name(d)) only accepts LinearCurve for loss models.", + ) + end + l1 = PSY.get_proportional_term(loss) + l0 = PSY.get_constant_term(loss) name = PSY.get_name(d) for t in get_time_steps(container) if l1 == 0.0 && l0 == 0.0 From 238d9b73b726e20245cb4ca5b77a385d0cb6d781 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Thu, 8 Aug 2024 20:50:18 -0700 Subject: [PATCH 02/68] add PWL Flow Calculations --- .../devices/TwoTerminalDC_branches.jl | 122 ++++++++++++++++++ 1 file changed, 122 insertions(+) diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 97c693635d..7ea74b76df 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -151,6 +151,128 @@ get_initial_conditions_device_model( ) where {T <: PSY.TwoTerminalHVDCLine, U <: AbstractTwoTerminalDCLineFormulation} = DeviceModel(T, U) +####################################### PWL Constraints ####################################################### + +function _get_pwl_loss_params(d::PSY.TwoTerminalHVDCLine, loss::PSY.LinearCurve) + from_to_loss_params = Vector{Float64}(undef, 4) + to_from_loss_params = Vector{Float64}(undef, 4) + loss_factor = PSY.get_proportional_term(loss) + P_send0 = PSY.get_constant_term(loss) + P_max_ft = PSY.get_active_power_limits_from(d).max + P_max_tf = PSY.get_active_power_limits_to(d).max + if P_max_ft != P_max_tf + error( + "HVDC Line $(PSY.get_name(d)) has non-symmetrical limits for from and to, that are not supported in the HVDCTwoTerminalPiecewiseLoss formulation", + ) + end + P_sendS = P_max_ft + ### Update Params Vectors ### + from_to_loss_params[1] = -P_sendS - P_send0 + from_to_loss_params[2] = -P_send0 + from_to_loss_params[3] = 0.0 + from_to_loss_params[4] = P_sendS * (1 - loss_factor) + + to_from_loss_params[1] = P_sendS * (1 - loss_factor) + to_from_loss_params[2] = 0.0 + to_from_loss_params[3] = -P_send0 + to_from_loss_params[4] = -P_sendS - P_send0 + + return from_to_loss_params, to_from_loss_params +end + +function _get_pwl_loss_params( + d::PSY.TwoTerminalHVDCLine, + loss::PSY.PiecewiseIncrementalCurve, +) + p_breakpoints = PSY.get_x_coords(loss) + loss_factors = PSY.get_slopes(loss) + len_segments = length(loss_factors) + len_variables = 2 * len_segments + 2 + from_to_loss_params = Vector{Float64}(undef, len_variables) + to_from_loss_params = similar(from_to_loss_params) + P_max_ft = PSY.get_active_power_limits_from(d).max + P_max_tf = PSY.get_active_power_limits_to(d).max + if P_max_ft != P_max_tf + error( + "HVDC Line $(PSY.get_name(d)) has non-symmetrical limits for from and to, that are not supported in the HVDCTwoTerminalPiecewiseLoss formulation", + ) + end + if P_max_ft != last(p_breakpoints) + error( + "Maximum power limit $P_max_ft of HVDC Line $(PSY.get_name(d)) has different value of last breakpoint from Loss data $(last(p_breakpoints)).", + ) + end + ### Update Params Vectors ### + ## Update from 1 to S + for i in 1:len_segments + from_to_loss_params[i] = -p_breakpoints[2 + len_segments - i] - p_breakpoints[1] # for i = 1: P_end, for i = len_segments: P_2 + to_from_loss_params[i] = + p_breakpoints[2 + len_segments - i] * (1 - loss_factors[len_segments + 1 - i]) + end + ## Update from S+1 and S+2 + from_to_loss_params[len_segments + 1] = -p_breakpoints[1] # P_send0 + from_to_loss_params[len_segments + 2] = 0.0 + to_from_loss_params[len_segments + 1] = 0.0 + to_from_loss_params[len_segments + 2] = -p_breakpoints[1] # P_send0 + ## Update from S+3 to 2S+2 + for i in 1:len_segments + from_to_loss_params[2 + len_segments + i] = + p_breakpoints[i + 1] * (1 - loss_factors[i]) + to_from_loss_params[2 + len_segments + i] = -p_breakpoints[i + 1] - p_breakpoints[1] + end + + return from_to_loss_params, to_from_loss_params +end + +function add_constraints!( + container::OptimizationContainer, + ::Type{T}, + devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}}, + ::DeviceModel{U, HVDCTwoTerminalPiecewiseLoss}, + ::NetworkModel{<:PM.AbstractPowerModel}, +) where {T <: HVDCFlowCalculationConstraint, U <: PSY.TwoTerminalHVDCLine} + var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) + names = axes(var_pwl)[1] + segments = axes(var_pwl)[2] + time_steps = get_time_steps(container) + flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) + flow_tf = get_variable(container, FlowActivePowerToFromVariable(), U) + + constraint_from_to = + add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") + constraint_to_from = + add_constraints_container!(container, T(), U, names, time_steps; meta = "tf") + for d in devices + name = PSY.get_name(d) + loss = PSY.get_loss(d) + from_to_params, to_from_params = _get_pwl_loss_params(d, loss) + for t in time_steps + ## Add Equality Constraints ## + constraint_from_to[PSY.get_name(d), t] = JuMP.@constraint( + get_jump_model(container), + flow_ft[name, t] == sum( + var_pwl[name, s, t] * from_to_params[ix] for + (ix, s) in enumerate(segments) + ) + ) + constraint_to_from[PSY.get_name(d), t] = JuMP.@constraint( + get_jump_model(container), + flow_tf[name, t] == sum( + var_pwl[name, s, t] * to_from_params[ix] for + (ix, s) in enumerate(segments) + ) + ) + ## Add SOS Constraints ### + pwl_vars_subset = [var_pwl[name, s, t] for s in segments] + JuMP.@constraint( + get_jump_model(container), + pwl_vars_subset in MOI.SOS2(collect(1:length(segments))) + ) + end + end + return +end + #################################### Rate Limits Constraints ################################################## function _get_flow_bounds(d::PSY.TwoTerminalHVDCLine) check_hvdc_line_limits_consistency(d) From 495a168700e70f0556cbc3e819aac0c66b55bff0 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Thu, 8 Aug 2024 20:52:27 -0700 Subject: [PATCH 03/68] comment TODO --- src/devices_models/device_constructors/branch_constructor.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 04767cb3f4..b1585d4e3c 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -819,6 +819,8 @@ function construct_device!( return end +# TODO: Other models besides PTDF +#= function construct_device!( container::OptimizationContainer, sys::PSY.System, @@ -876,6 +878,7 @@ function construct_device!( add_constraint_dual!(container, sys, model) return end +=# function construct_device!( container::OptimizationContainer, From a76e6400e4c5f9d4bff5aba3d3c1221908bfa3df Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Thu, 8 Aug 2024 20:56:58 -0700 Subject: [PATCH 04/68] add LL and UB for SOS variable --- src/devices_models/devices/TwoTerminalDC_branches.jl | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 7ea74b76df..55e75517a7 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -131,6 +131,18 @@ function get_variable_upper_bound( end end +get_variable_upper_bound( + ::HVDCPiecewiseLossVariable, + d::PSY.TwoTerminalHVDCLine, + ::Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, +) = 1.0 + +get_variable_lower_bound( + ::HVDCPiecewiseLossVariable, + d::PSY.TwoTerminalHVDCLine, + ::Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, +) = 0.0 + function get_default_time_series_names( ::Type{U}, ::Type{V}, From 7e08060988617ca74dbcbdd76d4f38fb50aef326 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 13 Aug 2024 20:08:19 -0700 Subject: [PATCH 05/68] add sparse method for PWL vars --- src/core/variables.jl | 7 +-- .../device_constructors/branch_constructor.jl | 3 +- src/devices_models/devices/AC_branches.jl | 47 +++++++++++++++++++ .../devices/TwoTerminalDC_branches.jl | 23 ++++++--- 4 files changed, 69 insertions(+), 11 deletions(-) diff --git a/src/core/variables.jl b/src/core/variables.jl index 53fe2de1b8..4b746b26ae 100644 --- a/src/core/variables.jl +++ b/src/core/variables.jl @@ -228,14 +228,14 @@ Docs abbreviation: ``u^\\text{dir}`` """ struct HVDCFlowDirectionVariable <: VariableType end +abstract type SparseVariableType <: VariableType end + """ Struct to dispatch the creation of HVDC Piecewise Loss Variables Docs abbreviation: ``h`` """ -struct HVDCPiecewiseLossVariable <: VariableType end - -abstract type SparseVariableType <: VariableType end +struct HVDCPiecewiseLossVariable <: SparseVariableType end """ Struct to dispatch the creation of piecewise linear cost variables for objective function @@ -281,6 +281,7 @@ const START_VARIABLES = (HotStartVariable, WarmStartVariable, ColdStartVariable) should_write_resulting_value(::Type{PieceWiseLinearCostVariable}) = false should_write_resulting_value(::Type{PieceWiseLinearBlockOffer}) = false +should_write_resulting_value(::Type{HVDCPiecewiseLossVariable}) = false convert_result_to_natural_units(::Type{ActivePowerVariable}) = true convert_result_to_natural_units(::Type{PowerAboveMinimumVariable}) = true convert_result_to_natural_units(::Type{ActivePowerInVariable}) = true diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index b1585d4e3c..e5349fcc4d 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -778,7 +778,8 @@ function construct_device!( devices, HVDCTwoTerminalPiecewiseLoss(), ) - _add_dense_pwl_loss_variables!(container, devices, model) + #_add_dense_pwl_loss_variables!(container, devices, model) + _add_sparse_pwl_loss_variables!(container, devices, model) add_to_expression!( container, ActivePowerBalance, diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 06eb06536c..93e7e15fab 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -248,6 +248,53 @@ function _add_dense_pwl_loss_variables!( end end +function _add_sparse_pwl_loss_variables!( + container::OptimizationContainer, + devices, + model::DeviceModel{D, HVDCTwoTerminalPiecewiseLoss}, +) where {D <: TwoTerminalHVDCTypes} + # Check if type and length of PWL loss model are the same for all devices + #_check_pwl_loss_model(devices) + + # Create Variables + time_steps = get_time_steps(container) + settings = get_settings(container) + formulation = HVDCTwoTerminalPiecewiseLoss() + T = HVDCPiecewiseLossVariable + binary = get_variable_binary(T(), D, formulation) + first_loss = PSY.get_loss(first(devices)) + if isa(first_loss, PSY.LinearCurve) + len_segments = 4 # 2*1 + 2 + elseif isa(first_loss, PSY.PiecewiseIncrementalCurve) + len_segments = 2 * length(PSY.get_slopes(first_loss)) + 2 + else + error("Should not be here") + end + + T = HVDCPiecewiseLossVariable + var_container = lazy_container_addition!(container, T(), D) + + for d in devices + name = PSY.get_name(d) + for t in time_steps + pwlvars = Array{JuMP.VariableRef}(undef, len_segments) + for i in 1:len_segments + pwlvars[i] = + var_container[(name, i, t)] = JuMP.@variable( + get_jump_model(container), + base_name = "$(T)_$(name)_{pwl_$(i), $(t)}", + binary = binary + ) + ub = get_variable_upper_bound(T(), d, formulation) + ub !== nothing && JuMP.set_upper_bound(var_container[name, i, t], ub) + + lb = get_variable_lower_bound(T(), d, formulation) + lb !== nothing && JuMP.set_lower_bound(var_container[name, i, t], lb) + end + end + end +end + ################################## Rate Limits constraint_infos ############################ """ diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 55e75517a7..236eb005df 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -165,6 +165,15 @@ get_initial_conditions_device_model( ####################################### PWL Constraints ####################################################### +function _get_range_segments(::PSY.TwoTerminalHVDCLine, loss::PSY.LinearCurve) + return 1:4 +end + +function _get_range_segments(::PSY.TwoTerminalHVDCLine, loss::PSY.PiecewiseIncrementalCurve) + loss_factors = PSY.get_slopes(loss) + return 1:length(loss_factors) +end + function _get_pwl_loss_params(d::PSY.TwoTerminalHVDCLine, loss::PSY.LinearCurve) from_to_loss_params = Vector{Float64}(undef, 4) to_from_loss_params = Vector{Float64}(undef, 4) @@ -244,8 +253,7 @@ function add_constraints!( ::NetworkModel{<:PM.AbstractPowerModel}, ) where {T <: HVDCFlowCalculationConstraint, U <: PSY.TwoTerminalHVDCLine} var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) - names = axes(var_pwl)[1] - segments = axes(var_pwl)[2] + names = PSY.get_name.(devices) time_steps = get_time_steps(container) flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) flow_tf = get_variable(container, FlowActivePowerToFromVariable(), U) @@ -258,27 +266,28 @@ function add_constraints!( name = PSY.get_name(d) loss = PSY.get_loss(d) from_to_params, to_from_params = _get_pwl_loss_params(d, loss) + segments = _get_range_segments(d, loss) for t in time_steps ## Add Equality Constraints ## constraint_from_to[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), flow_ft[name, t] == sum( - var_pwl[name, s, t] * from_to_params[ix] for - (ix, s) in enumerate(segments) + var_pwl[name, ix, t] * from_to_params[ix] for + ix in segments ) ) constraint_to_from[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), flow_tf[name, t] == sum( - var_pwl[name, s, t] * to_from_params[ix] for - (ix, s) in enumerate(segments) + var_pwl[name, ix, t] * to_from_params[ix] for + ix in segments ) ) ## Add SOS Constraints ### pwl_vars_subset = [var_pwl[name, s, t] for s in segments] JuMP.@constraint( get_jump_model(container), - pwl_vars_subset in MOI.SOS2(collect(1:length(segments))) + pwl_vars_subset in MOI.SOS2(collect(segments)) ) end end From aa97451f87bc4e73b4fce38f6b83ccd13e98457a Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 14 Aug 2024 19:55:47 -0700 Subject: [PATCH 06/68] add SOS formulation --- src/core/formulations.jl | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/core/formulations.jl b/src/core/formulations.jl index 5408d934be..49bd495d71 100644 --- a/src/core/formulations.jl +++ b/src/core/formulations.jl @@ -138,6 +138,11 @@ struct HVDCTwoTerminalDispatch <: AbstractTwoTerminalDCLineFormulation end Branch type to represent piecewise lossy power flow on two terminal DC lines """ struct HVDCTwoTerminalPiecewiseLoss <: AbstractTwoTerminalDCLineFormulation end +""" +Branch type to represent piecewise lossy power flow on two terminal DC lines using SOS2 +""" +struct HVDCTwoTerminalSOSPiecewiseLoss <: AbstractTwoTerminalDCLineFormulation end + # Not Implemented # struct VoltageSourceDC <: AbstractTwoTerminalDCLineFormulation end From d46bd2119eb0f413a731643a2cf13aad374e28e3 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 14 Aug 2024 19:55:56 -0700 Subject: [PATCH 07/68] fix add to expression for PWL model --- .../devices/common/add_to_expression.jl | 37 +++++++++++++++++++ 1 file changed, 37 insertions(+) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index eb7702ed08..d839a40f28 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -375,6 +375,43 @@ function add_to_expression!( return end +""" +PWL implementation to add FromTo branch variables to SystemBalanceExpressions +""" +function add_to_expression!( + container::OptimizationContainer, + ::Type{T}, + ::Type{U}, + devices::IS.FlattenIteratorWrapper{V}, + ::DeviceModel{V, W}, + network_model::NetworkModel{X}, +) where { + T <: ActivePowerBalance, + U <: FlowActivePowerFromToVariable, + V <: TwoTerminalHVDCTypes, + W <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + X <: AbstractPTDFModel, +} + var = get_variable(container, U(), V) + nodal_expr = get_expression(container, T(), PSY.ACBus) + sys_expr = get_expression(container, T(), _system_expression_type(X)) + 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) + 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) + flow_variable = var[PSY.get_name(d), t] + _add_to_jump_expression!(nodal_expr[bus_no_from, t], flow_variable, 1.0) + if ref_bus_from != ref_bus_to + _add_to_jump_expression!(sys_expr[ref_bus_from, t], flow_variable, 1.0) + end + end + end + return +end + """ Default implementation to add branch variables to SystemBalanceExpressions """ From 9203617589d62060b8c0e609cf9acf2d3b8df0d3 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 14 Aug 2024 19:56:19 -0700 Subject: [PATCH 08/68] add Z,W model --- src/core/variables.jl | 10 ++- .../device_constructors/branch_constructor.jl | 14 ++- src/devices_models/devices/AC_branches.jl | 62 ++++++++++++- .../devices/TwoTerminalDC_branches.jl | 88 ++++++++++++++++++- 4 files changed, 165 insertions(+), 9 deletions(-) diff --git a/src/core/variables.jl b/src/core/variables.jl index 4b746b26ae..75b927c36d 100644 --- a/src/core/variables.jl +++ b/src/core/variables.jl @@ -233,10 +233,17 @@ abstract type SparseVariableType <: VariableType end """ Struct to dispatch the creation of HVDC Piecewise Loss Variables -Docs abbreviation: ``h`` +Docs abbreviation: ``h`` or ``w`` """ struct HVDCPiecewiseLossVariable <: SparseVariableType end +""" +Struct to dispatch the creation of HVDC Piecewise Binary Loss Variables + +Docs abbreviation: ``z`` +""" +struct HVDCPiecewiseBinaryLossVariable <: SparseVariableType end + """ Struct to dispatch the creation of piecewise linear cost variables for objective function @@ -282,6 +289,7 @@ const START_VARIABLES = (HotStartVariable, WarmStartVariable, ColdStartVariable) should_write_resulting_value(::Type{PieceWiseLinearCostVariable}) = false should_write_resulting_value(::Type{PieceWiseLinearBlockOffer}) = false should_write_resulting_value(::Type{HVDCPiecewiseLossVariable}) = false +should_write_resulting_value(::Type{HVDCPiecewiseBinaryLossVariable}) = false convert_result_to_natural_units(::Type{ActivePowerVariable}) = true convert_result_to_natural_units(::Type{PowerAboveMinimumVariable}) = true convert_result_to_natural_units(::Type{ActivePowerInVariable}) = true diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index e5349fcc4d..5b1e232692 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -761,9 +761,12 @@ function construct_device!( container::OptimizationContainer, sys::PSY.System, ::ArgumentConstructStage, - model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + model::DeviceModel{T, U}, network_model::NetworkModel{<:AbstractPTDFModel}, -) where {T <: TwoTerminalHVDCTypes} +) where { + T <: TwoTerminalHVDCTypes, + U <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, +} devices = get_available_components(model, sys) add_variables!( @@ -803,9 +806,12 @@ function construct_device!( container::OptimizationContainer, sys::PSY.System, ::ModelConstructStage, - model::DeviceModel{T, HVDCTwoTerminalPiecewiseLoss}, + model::DeviceModel{T, U}, network_model::NetworkModel{<:AbstractPTDFModel}, -) where {T <: TwoTerminalHVDCTypes} +) where { + T <: TwoTerminalHVDCTypes, + U <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, +} devices = get_available_components(model, sys) add_constraints!(container, FlowRateConstraintFromTo, devices, model, network_model) diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 93e7e15fab..258ed0248c 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -33,6 +33,7 @@ get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.ACBranch}, :: get_variable_binary(::FlowActivePowerSlackUpperBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::HVDCPiecewiseLossVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false +get_variable_binary(::HVDCPiecewiseBinaryLossVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = true get_variable_upper_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = nothing get_variable_lower_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = 0.0 get_variable_upper_bound(::FlowActivePowerSlackLowerBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = nothing @@ -248,10 +249,11 @@ function _add_dense_pwl_loss_variables!( end end +#SOS Model function _add_sparse_pwl_loss_variables!( container::OptimizationContainer, devices, - model::DeviceModel{D, HVDCTwoTerminalPiecewiseLoss}, + ::DeviceModel{D, HVDCTwoTerminalSOSPiecewiseLoss}, ) where {D <: TwoTerminalHVDCTypes} # Check if type and length of PWL loss model are the same for all devices #_check_pwl_loss_model(devices) @@ -295,6 +297,64 @@ function _add_sparse_pwl_loss_variables!( end end +# Full Binary +function _add_sparse_pwl_loss_variables!( + container::OptimizationContainer, + devices, + ::DeviceModel{D, HVDCTwoTerminalPiecewiseLoss}, +) where {D <: TwoTerminalHVDCTypes} + # Check if type and length of PWL loss model are the same for all devices + #_check_pwl_loss_model(devices) + + # Create Variables + time_steps = get_time_steps(container) + settings = get_settings(container) + formulation = HVDCTwoTerminalPiecewiseLoss() + T = HVDCPiecewiseLossVariable + binary_T = get_variable_binary(T(), D, formulation) + U = HVDCPiecewiseBinaryLossVariable + binary_U = get_variable_binary(U(), D, formulation) + first_loss = PSY.get_loss(first(devices)) + if isa(first_loss, PSY.LinearCurve) + len_segments = 3 # 2*1 + 1 + elseif isa(first_loss, PSY.PiecewiseIncrementalCurve) + len_segments = 2 * length(PSY.get_slopes(first_loss)) + 1 + else + error("Should not be here") + end + + var_container = lazy_container_addition!(container, T(), D) + var_container_binary = lazy_container_addition!(container, U(), D) + + for d in devices + name = PSY.get_name(d) + for t in time_steps + pwlvars = Array{JuMP.VariableRef}(undef, len_segments) + pwlvars_bin = Array{JuMP.VariableRef}(undef, len_segments) + for i in 1:len_segments + pwlvars[i] = + var_container[(name, i, t)] = JuMP.@variable( + get_jump_model(container), + base_name = "$(T)_$(name)_{pwl_$(i), $(t)}", + binary = binary_T + ) + ub = get_variable_upper_bound(T(), d, formulation) + ub !== nothing && JuMP.set_upper_bound(var_container[name, i, t], ub) + + lb = get_variable_lower_bound(T(), d, formulation) + lb !== nothing && JuMP.set_lower_bound(var_container[name, i, t], lb) + + pwlvars_bin[i] = + var_container_binary[(name, i, t)] = JuMP.@variable( + get_jump_model(container), + base_name = "$(U)_$(name)_{pwl_$(i), $(t)}", + binary = binary_U + ) + end + end + end +end + ################################## Rate Limits constraint_infos ############################ """ diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 236eb005df..1bad48cfee 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -171,7 +171,7 @@ end function _get_range_segments(::PSY.TwoTerminalHVDCLine, loss::PSY.PiecewiseIncrementalCurve) loss_factors = PSY.get_slopes(loss) - return 1:length(loss_factors) + return 1:(2 * length(loss_factors) + 2) end function _get_pwl_loss_params(d::PSY.TwoTerminalHVDCLine, loss::PSY.LinearCurve) @@ -253,6 +253,7 @@ function add_constraints!( ::NetworkModel{<:PM.AbstractPowerModel}, ) where {T <: HVDCFlowCalculationConstraint, U <: PSY.TwoTerminalHVDCLine} var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) + var_pwl_bin = get_variable(container, HVDCPiecewiseBinaryLossVariable(), U) names = PSY.get_name.(devices) time_steps = get_time_steps(container) flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) @@ -262,10 +263,83 @@ function add_constraints!( add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") constraint_to_from = add_constraints_container!(container, T(), U, names, time_steps; meta = "tf") + constraint_binary = + add_constraints_container!(container, T(), U, names, time_steps; meta = "bin") for d in devices name = PSY.get_name(d) loss = PSY.get_loss(d) from_to_params, to_from_params = _get_pwl_loss_params(d, loss) + #@show from_to_params + #@show to_from_params + range_segments = 1:(length(from_to_params) - 1) # 1:(2S+1) + for t in time_steps + ## Add Equality Constraints ## + constraint_from_to[name, t] = JuMP.@constraint( + get_jump_model(container), + flow_ft[name, t] == + sum( + var_pwl_bin[name, ix, t] * from_to_params[ix] for + ix in range_segments + ) + sum( + var_pwl[name, ix, t] * (from_to_params[ix + 1] - from_to_params[ix]) for + ix in range_segments + ) + ) + constraint_to_from[name, t] = JuMP.@constraint( + get_jump_model(container), + flow_tf[name, t] == + sum( + var_pwl_bin[name, ix, t] * to_from_params[ix] for + ix in range_segments + ) + sum( + var_pwl[name, ix, t] * (to_from_params[ix + 1] - to_from_params[ix]) for + ix in range_segments + ) + ) + ## Add Binary Bound ### + constraint_binary[name, t] = JuMP.@constraint( + get_jump_model(container), + sum(var_pwl_bin[name, ix, t] for ix in range_segments) == 1.0 + ) + ## Add Bounds for Continuous ## + for ix in range_segments + JuMP.@constraint( + get_jump_model(container), + var_pwl[name, ix, t] <= var_pwl_bin[name, ix, t] + ) + if ix == div(length(range_segments) + 1, 2) + JuMP.fix(var_pwl[name, ix, t], 0.0; force = true) + end + end + end + end + return +end + +# SOS Model +function add_constraints!( + container::OptimizationContainer, + ::Type{T}, + devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}}, + ::DeviceModel{U, HVDCTwoTerminalSOSPiecewiseLoss}, + ::NetworkModel{<:PM.AbstractPowerModel}, +) where {T <: HVDCFlowCalculationConstraint, U <: PSY.TwoTerminalHVDCLine} + var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) + names = PSY.get_name.(devices) + time_steps = get_time_steps(container) + flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) + flow_tf = get_variable(container, FlowActivePowerToFromVariable(), U) + + constraint_from_to = + add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") + constraint_to_from = + add_constraints_container!(container, T(), U, names, time_steps; meta = "tf") + for d in devices + name = PSY.get_name(d) + loss = PSY.get_loss(d) + from_to_params, to_from_params = _get_pwl_loss_params(d, loss) + #@show from_to_params + #@show to_from_params segments = _get_range_segments(d, loss) for t in time_steps ## Add Equality Constraints ## @@ -439,7 +513,11 @@ function add_constraints!( ) where { T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, + V <: Union{ + HVDCTwoTerminalDispatch, + HVDCTwoTerminalPiecewiseLoss, + HVDCTwoTerminalSOSPiecewiseLoss, + }, } inter_network_branches = U[] for d in devices @@ -464,7 +542,11 @@ function add_constraints!( ) where { T <: Union{FlowRateConstraintToFrom, FlowRateConstraintFromTo}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{HVDCTwoTerminalDispatch, HVDCTwoTerminalPiecewiseLoss}, + V <: Union{ + HVDCTwoTerminalDispatch, + HVDCTwoTerminalPiecewiseLoss, + HVDCTwoTerminalSOSPiecewiseLoss, + }, } _add_hvdc_flow_constraints!(container, devices, T()) return From 9cc16bfabbe4e567e379b5a9d0fdede6db06655d Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 10 Sep 2024 09:25:31 -0700 Subject: [PATCH 09/68] remove compact stuff --- .../common/objective_function/market_bid.jl | 15 --------------- 1 file changed, 15 deletions(-) diff --git a/src/devices_models/devices/common/objective_function/market_bid.jl b/src/devices_models/devices/common/objective_function/market_bid.jl index d46025e62f..635f59e07d 100644 --- a/src/devices_models/devices/common/objective_function/market_bid.jl +++ b/src/devices_models/devices/common/objective_function/market_bid.jl @@ -302,21 +302,6 @@ function _add_pwl_term!( device_base_power, ) - compact_status = validate_compact_pwl_data(component, data, base_power) - if !uses_compact_power(component, V()) && compact_status == COMPACT_PWL_STATUS.VALID - error( - "The data provided is not compatible with formulation $V. Use a formulation compatible with Compact Cost Functions", - ) - # data = _convert_to_full_variable_cost(data, component) - elseif uses_compact_power(component, V()) && compact_status != COMPACT_PWL_STATUS.VALID - @warn( - "The cost data provided is not in compact form. Will attempt to convert. Errors may occur." - ) - data = convert_to_compact_variable_cost(data) - else - @debug uses_compact_power(component, V()) compact_status name T V - end - cost_is_convex = PSY.is_convex(data) if !cost_is_convex error("MarketBidCost for component $(name) is non-convex") From 477a98bd8e2a5e67276c11b5968d159906a23d2b Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 11 Sep 2024 15:48:09 -0700 Subject: [PATCH 10/68] add better variable names --- src/core/variables.jl | 14 ++ .../device_constructors/branch_constructor.jl | 9 +- src/devices_models/devices/AC_branches.jl | 2 + .../devices/TwoTerminalDC_branches.jl | 121 +++++++++++++++--- .../devices/common/add_to_expression.jl | 43 ++++++- 5 files changed, 164 insertions(+), 25 deletions(-) diff --git a/src/core/variables.jl b/src/core/variables.jl index 75b927c36d..3050c31f51 100644 --- a/src/core/variables.jl +++ b/src/core/variables.jl @@ -228,6 +228,20 @@ Docs abbreviation: ``u^\\text{dir}`` """ struct HVDCFlowDirectionVariable <: VariableType end +""" +Struct to dispatch the creation of HVDC Received Flow at From Bus Variables for PWL formulations + +Docs abbreviation: ``x`` +""" +struct HVDCActivePowerReceivedFromVariable <: VariableType end + +""" +Struct to dispatch the creation of HVDC Received Flow at To Bus Variables for PWL formulations + +Docs abbreviation: ``y`` +""" +struct HVDCActivePowerReceivedToVariable <: VariableType end + abstract type SparseVariableType <: VariableType end """ diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 5b1e232692..94c6f4bc61 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -771,22 +771,21 @@ function construct_device!( get_available_components(model, sys) add_variables!( container, - FlowActivePowerToFromVariable, + HVDCActivePowerReceivedFromVariable, devices, HVDCTwoTerminalPiecewiseLoss(), ) add_variables!( container, - FlowActivePowerFromToVariable, + HVDCActivePowerReceivedToVariable, devices, HVDCTwoTerminalPiecewiseLoss(), ) - #_add_dense_pwl_loss_variables!(container, devices, model) _add_sparse_pwl_loss_variables!(container, devices, model) add_to_expression!( container, ActivePowerBalance, - FlowActivePowerToFromVariable, + HVDCActivePowerReceivedFromVariable, devices, model, network_model, @@ -794,7 +793,7 @@ function construct_device!( add_to_expression!( container, ActivePowerBalance, - FlowActivePowerFromToVariable, + HVDCActivePowerReceivedToVariable, devices, model, network_model, diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index 258ed0248c..e20d1edbfc 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -33,6 +33,8 @@ get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.ACBranch}, :: get_variable_binary(::FlowActivePowerSlackUpperBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::FlowActivePowerSlackLowerBound, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::HVDCPiecewiseLossVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false +get_variable_binary(::HVDCActivePowerReceivedFromVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false +get_variable_binary(::HVDCActivePowerReceivedToVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = false get_variable_binary(::HVDCPiecewiseBinaryLossVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation,) = true get_variable_upper_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = nothing get_variable_lower_bound(::FlowActivePowerSlackUpperBound, ::PSY.ACBranch, ::AbstractBranchFormulation) = 0.0 diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 1bad48cfee..26acd94857 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -111,6 +111,30 @@ get_variable_lower_bound( ::HVDCTwoTerminalDispatch, ) = PSY.get_active_power_limits_to(d).min +get_variable_upper_bound( + ::HVDCActivePowerReceivedFromVariable, + d::PSY.TwoTerminalHVDCLine, + ::AbstractTwoTerminalDCLineFormulation, +) = PSY.get_active_power_limits_from(d).max + +get_variable_lower_bound( + ::HVDCActivePowerReceivedFromVariable, + d::PSY.TwoTerminalHVDCLine, + ::AbstractTwoTerminalDCLineFormulation, +) = PSY.get_active_power_limits_from(d).min + +get_variable_upper_bound( + ::HVDCActivePowerReceivedToVariable, + d::PSY.TwoTerminalHVDCLine, + ::AbstractTwoTerminalDCLineFormulation, +) = PSY.get_active_power_limits_to(d).max + +get_variable_lower_bound( + ::HVDCActivePowerReceivedToVariable, + d::PSY.TwoTerminalHVDCLine, + ::AbstractTwoTerminalDCLineFormulation, +) = PSY.get_active_power_limits_to(d).min + function get_variable_upper_bound( ::HVDCLosses, d::PSY.TwoTerminalHVDCLine, @@ -256,8 +280,8 @@ function add_constraints!( var_pwl_bin = get_variable(container, HVDCPiecewiseBinaryLossVariable(), U) names = PSY.get_name.(devices) time_steps = get_time_steps(container) - flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) - flow_tf = get_variable(container, FlowActivePowerToFromVariable(), U) + flow_ft = get_variable(container, HVDCActivePowerReceivedFromVariable(), U) + flow_tf = get_variable(container, HVDCActivePowerReceivedToVariable(), U) constraint_from_to = add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") @@ -327,8 +351,8 @@ function add_constraints!( var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) names = PSY.get_name.(devices) time_steps = get_time_steps(container) - flow_ft = get_variable(container, FlowActivePowerFromToVariable(), U) - flow_tf = get_variable(container, FlowActivePowerToFromVariable(), U) + flow_ft = get_variable(container, HVDCActivePowerReceivedFromVariable(), U) + flow_tf = get_variable(container, HVDCActivePowerReceivedToVariable(), U) constraint_from_to = add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") @@ -475,7 +499,12 @@ end function _add_hvdc_flow_constraints!( container::OptimizationContainer, devices::Union{Vector{T}, IS.FlattenIteratorWrapper{T}}, - var::Union{FlowActivePowerFromToVariable, FlowActivePowerToFromVariable}, + var::Union{ + FlowActivePowerFromToVariable, + FlowActivePowerToFromVariable, + HVDCActivePowerReceivedFromVariable, + HVDCActivePowerReceivedToVariable, + }, constraint::Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, ) where {T <: PSY.TwoTerminalHVDCLine} time_steps = get_time_steps(container) @@ -508,16 +537,11 @@ function add_constraints!( container::OptimizationContainer, ::Type{T}, devices::IS.FlattenIteratorWrapper{U}, - model::DeviceModel{U, V}, + model::DeviceModel{U, HVDCTwoTerminalDispatch}, network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{ - HVDCTwoTerminalDispatch, - HVDCTwoTerminalPiecewiseLoss, - HVDCTwoTerminalSOSPiecewiseLoss, - }, } inter_network_branches = U[] for d in devices @@ -537,16 +561,11 @@ function add_constraints!( container::OptimizationContainer, ::Type{T}, devices::IS.FlattenIteratorWrapper{U}, - ::DeviceModel{U, V}, + ::DeviceModel{U, HVDCTwoTerminalDispatch}, ::NetworkModel{<:PM.AbstractDCPModel}, ) where { T <: Union{FlowRateConstraintToFrom, FlowRateConstraintFromTo}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{ - HVDCTwoTerminalDispatch, - HVDCTwoTerminalPiecewiseLoss, - HVDCTwoTerminalSOSPiecewiseLoss, - }, } _add_hvdc_flow_constraints!(container, devices, T()) return @@ -566,6 +585,74 @@ function add_constraints!( return end +function add_constraints!( + container::OptimizationContainer, + ::Type{T}, + devices::IS.FlattenIteratorWrapper{U}, + model::DeviceModel{U, V}, + network_model::NetworkModel{CopperPlatePowerModel}, +) where { + T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, + U <: PSY.TwoTerminalHVDCLine, + V <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, +} + inter_network_branches = U[] + for d in devices + 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) + if ref_bus_from != ref_bus_to + push!(inter_network_branches, d) + end + end + if !isempty(inter_network_branches) + if T <: FlowRateConstraintFromTo + _add_hvdc_flow_constraints!( + container, + devices, + HVDCActivePowerReceivedFromVariable(), + T(), + ) + else + _add_hvdc_flow_constraints!( + container, + devices, + HVDCActivePowerReceivedToVariable(), + T(), + ) + end + end + return +end + +function add_constraints!( + container::OptimizationContainer, + ::Type{T}, + devices::IS.FlattenIteratorWrapper{U}, + ::DeviceModel{U, V}, + ::NetworkModel{<:AbstractPTDFModel}, +) where { + T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, + U <: PSY.TwoTerminalHVDCLine, + V <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, +} + if T <: FlowRateConstraintFromTo + _add_hvdc_flow_constraints!( + container, + devices, + HVDCActivePowerReceivedFromVariable(), + T(), + ) + else + _add_hvdc_flow_constraints!( + container, + devices, + HVDCActivePowerReceivedToVariable(), + T(), + ) + end + return +end + function add_constraints!( container::OptimizationContainer, ::Type{T}, diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index d839a40f28..48d77cc8af 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -366,9 +366,9 @@ function add_to_expression!( ref_bus_from = get_reference_bus(network_model, PSY.get_arc(d).from) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] - _add_to_jump_expression!(nodal_expr[bus_no_from, t], flow_variable, -1.0) + _add_to_jump_expression!(nodal_expr[bus_no_from, t], flow_variable, 1.0) if ref_bus_from != ref_bus_to - _add_to_jump_expression!(sys_expr[ref_bus_from, t], flow_variable, -1.0) + _add_to_jump_expression!(sys_expr[ref_bus_from, t], flow_variable, 1.0) end end end @@ -387,7 +387,7 @@ function add_to_expression!( network_model::NetworkModel{X}, ) where { T <: ActivePowerBalance, - U <: FlowActivePowerFromToVariable, + U <: HVDCActivePowerReceivedFromVariable, V <: TwoTerminalHVDCTypes, W <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, X <: AbstractPTDFModel, @@ -412,6 +412,43 @@ function add_to_expression!( return end +""" +PWL implementation to add FromTo branch variables to SystemBalanceExpressions +""" +function add_to_expression!( + container::OptimizationContainer, + ::Type{T}, + ::Type{U}, + devices::IS.FlattenIteratorWrapper{V}, + ::DeviceModel{V, W}, + network_model::NetworkModel{X}, +) where { + T <: ActivePowerBalance, + U <: HVDCActivePowerReceivedToVariable, + V <: TwoTerminalHVDCTypes, + W <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + X <: AbstractPTDFModel, +} + var = get_variable(container, U(), V) + nodal_expr = get_expression(container, T(), PSY.ACBus) + sys_expr = get_expression(container, T(), _system_expression_type(X)) + 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_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) + flow_variable = var[PSY.get_name(d), t] + _add_to_jump_expression!(nodal_expr[bus_no_to, t], flow_variable, 1.0) + if ref_bus_from != ref_bus_to + _add_to_jump_expression!(sys_expr[ref_bus_to, t], flow_variable, 1.0) + end + end + end + return +end + """ Default implementation to add branch variables to SystemBalanceExpressions """ From 377cf9c357449b7b876fabab4f08f7ded30bdb2b Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 17 Sep 2024 14:21:39 -0600 Subject: [PATCH 11/68] fix addition --- src/core/optimization_container.jl | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index e3e7f2ffb7..10b7fa3eb8 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -33,7 +33,12 @@ function get_objective_expression(v::ObjectiveFunction) if iszero(v.variant_terms) return v.invariant_terms else - return JuMP.add_to_expression!(v.variant_terms, v.invariant_terms) + # JuMP doesn't support expression conversion from Affn to QuadExpressions + if isa(v.invariant_terms, JuMP.GenericQuadExpr) + return JuMP.add_to_expression!(v.invariant_terms, v.variant_terms) + else + return JuMP.add_to_expression!(v.variant_terms, v.invariant_terms) + end end end get_sense(v::ObjectiveFunction) = v.sense From 298a709947f36f214521f96289a32f57c647ae27 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 17 Sep 2024 13:46:13 -0700 Subject: [PATCH 12/68] add nonspin slack --- src/services_models/service_slacks.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/services_models/service_slacks.jl b/src/services_models/service_slacks.jl index 0ec2464fda..9d530e3e40 100644 --- a/src/services_models/service_slacks.jl +++ b/src/services_models/service_slacks.jl @@ -1,7 +1,7 @@ function reserve_slacks!( container::OptimizationContainer, service::T, -) where {T <: PSY.Reserve} +) where {T <: Union{PSY.Reserve, PSY.ReserveNonSpinning}} time_steps = get_time_steps(container) variable = add_variable_container!( container, From 9b20947a8d10652f04d3ef1e1ddda87070630555 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:15 -0600 Subject: [PATCH 13/68] save fixed value for variable --- src/core/optimization_container.jl | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index e3e7f2ffb7..c0fc9dfba6 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -1716,6 +1716,9 @@ function _process_duals(container::OptimizationContainer, lp_optimizer) if JuMP.has_upper_bound(first(variable)) cache[key][:ub] = JuMP.upper_bound.(variable) end + if JuMP.is_fixed(first(variable)) && is_integer_flag + cache[key][:fixed_int_value] = jump_value.(v) + end cache[key][:integer] = is_integer_flag JuMP.fix.(variable, var_cache[key]; force = true) end @@ -1756,6 +1759,9 @@ function _process_duals(container::OptimizationContainer, lp_optimizer) else JuMP.unfix.(variable) JuMP.set_binary.(variable) + if haskey(cache[key], :fixed_int_value) + JuMP.fix.(variable, cache[key][:fixed_int_value]) + end #= Needed if a model has integer variables if haskey(cache[key], :lb) && JuMP.has_lower_bound(first(variable)) JuMP.set_lower_bound.(variable, cache[key][:lb]) From 7d46f731e103c2ac02df653098fc478865d01495 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:34 -0600 Subject: [PATCH 14/68] change implementation of must_run vars --- .../devices/thermal_generation.jl | 28 ++++++++++--------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 28b12488da..709593521e 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -319,23 +319,25 @@ function add_variable!( container, variable_type, D, - [PSY.get_name(d) for d in devices], + [PSY.get_name(d) for d in devices if !PSY.get_must_run(d)], time_steps, ) - for t in time_steps, d in devices - name = PSY.get_name(d) - variable[name, t] = JuMP.@variable( - get_jump_model(container), - base_name = "$(T)_$(D)_{$(name), $(t)}", - binary = binary - ) - if get_warm_start(settings) - init = get_variable_warm_start_value(variable_type, d, formulation) - init !== nothing && JuMP.set_start_value(variable[name, t], init) - end + for d in devices if PSY.get_must_run(d) - JuMP.fix(variable[name, t], 1.0; force = true) + continue + end + name = PSY.get_name(d) + for t in time_steps + variable[name, t] = JuMP.@variable( + get_jump_model(container), + base_name = "$(T)_$(D)_{$(name), $(t)}", + binary = binary + ) + if get_warm_start(settings) + init = get_variable_warm_start_value(variable_type, d, formulation) + init !== nothing && JuMP.set_start_value(variable[name, t], init) + end end end From ba0bc69e7cfe6317197001d2ab0619434ff6b080 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:52 -0600 Subject: [PATCH 15/68] change implementation of range constraints with must run --- .../devices/common/range_constraint.jl | 102 ++++++++++++++---- 1 file changed, 82 insertions(+), 20 deletions(-) diff --git a/src/devices_models/devices/common/range_constraint.jl b/src/devices_models/devices/common/range_constraint.jl index a4e0566066..d92c5620e1 100644 --- a/src/devices_models/devices/common/range_constraint.jl +++ b/src/devices_models/devices/common/range_constraint.jl @@ -104,7 +104,7 @@ function _add_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_lb[ci_name, t] = - JuMP.@constraint(container.JuMPmodel, array[ci_name, t] >= limits.min) + JuMP.@constraint(get_jump_model(container), array[ci_name, t] >= limits.min) end return end @@ -126,7 +126,7 @@ function _add_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_ub[ci_name, t] = - JuMP.@constraint(container.JuMPmodel, array[ci_name, t] <= limits.max) + JuMP.@constraint(get_jump_model(container), array[ci_name, t] <= limits.max) end return end @@ -246,20 +246,86 @@ function _add_semicontinuous_lower_bound_range_constraints_impl!( ) where {T <: ConstraintType, V <: PSY.Component, W <: AbstractDeviceFormulation} time_steps = get_time_steps(container) names = [PSY.get_name(d) for d in devices] - binary_variables = [OnVariable()] + con_lb = add_constraints_container!(container, T(), V, names, time_steps; meta = "lb") + varbin = get_variable(container, OnVariable(), V) + + for device in devices + ci_name = PSY.get_name(device) + limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min * varbin[ci_name, t] + ) + end + end + return +end +function _add_semicontinuous_lower_bound_range_constraints_impl!( + container::OptimizationContainer, + ::Type{T}, + array, + devices::IS.FlattenIteratorWrapper{V}, + ::DeviceModel{V, W}, +) where {T <: ConstraintType, V <: PSY.ThermalGen, W <: AbstractDeviceFormulation} + time_steps = get_time_steps(container) + names = [PSY.get_name(d) for d in devices] con_lb = add_constraints_container!(container, T(), V, names, time_steps; meta = "lb") + varbin = get_variable(container, OnVariable(), V) - @assert length(binary_variables) == 1 "Expected $(binary_variables) for $U $V $T $W to be length 1" - varbin = get_variable(container, only(binary_variables), V) + for device in devices + ci_name = PSY.get_name(device) + limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type + if PSY.get_must_run_device(device) + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min + ) + end + else + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min * varbin[ci_name, t] + ) + end + end + end + return +end - for device in devices, t in time_steps +function _add_semicontinuous_upper_bound_range_constraints_impl!( + container::OptimizationContainer, + ::Type{T}, + array, + devices::IS.FlattenIteratorWrapper{V}, + model::DeviceModel{V, W}, +) where {T <: ConstraintType, V <: PSY.ThermalGen, W <: AbstractDeviceFormulation} + time_steps = get_time_steps(container) + names = [PSY.get_name(d) for d in devices] + con_ub = add_constraints_container!(container, T(), V, names, time_steps; meta = "ub") + varbin = get_variable(container, OnVariable(), V) + + for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, - array[ci_name, t] >= limits.min * varbin[ci_name, t] - ) + if PSY.get_must_run_device(device) + for t in time_steps + con_ub[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] <= limits.max + ) + end + else + for t in time_steps + con_ub[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] <= limits.max * varbin[ci_name, t] + ) + end + end end return end @@ -273,18 +339,14 @@ function _add_semicontinuous_upper_bound_range_constraints_impl!( ) where {T <: ConstraintType, V <: PSY.Component, W <: AbstractDeviceFormulation} time_steps = get_time_steps(container) names = [PSY.get_name(d) for d in devices] - binary_variables = [OnVariable()] - con_ub = add_constraints_container!(container, T(), V, names, time_steps; meta = "ub") - - @assert length(binary_variables) == 1 "Expected $(binary_variables) for $U $V $T $W to be length 1" - varbin = get_variable(container, only(binary_variables), V) + varbin = get_variable(container, OnVariable(), V) for device in devices, t in time_steps ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * varbin[ci_name, t] ) end @@ -375,7 +437,7 @@ function _add_reserve_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] >= limits.min * (1 - varbin[ci_name, t]) ) end @@ -409,7 +471,7 @@ function _add_reserve_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * (1 - varbin[ci_name, t]) ) end @@ -513,7 +575,7 @@ function _add_reserve_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, X) # depends on constraint type and formulation type con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] >= limits.min * varbin[ci_name, t] ) end @@ -545,7 +607,7 @@ function _add_reserve_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, X) # depends on constraint type and formulation type con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * varbin[ci_name, t] ) end From cc436e1a221a61914d15410e70cb93c6fe444d57 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:57:14 -0600 Subject: [PATCH 16/68] WP fix addition to expressions --- .../devices/common/add_to_expression.jl | 46 +++++++++++++------ 1 file changed, 32 insertions(+), 14 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index efbda4845c..3c7ce3277b 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -530,15 +530,24 @@ 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) - for d in devices, t in get_time_steps(container) + for d in devices name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) - _add_to_jump_expression!( - expression[bus_no, t], - variable[name, t], - get_variable_multiplier(U(), d, W()), - ) + for t in get_time_steps(container) + if PSY.get_must_run(d) + _add_to_jump_expression!( + expression[bus_no, t], + get_variable_multiplier(U(), d, W()), + ) + else + _add_to_jump_expression!( + expression[bus_no, t], + variable[name, t], + get_variable_multiplier(U(), d, W()), + ) + end + end end return end @@ -558,15 +567,24 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - for d in devices, t in get_time_steps(container) + for d in devices + name = PSY.get_name(d) bus = PSY.get_bus(d) area_name = PSY.get_name(PSY.get_area(bus)) - name = PSY.get_name(d) - _add_to_jump_expression!( - expression[area_name, t], - variable[name, t], - get_variable_multiplier(U(), d, W()), - ) + for t in get_time_steps(container) + if PSY.get_must_run(d) + _add_to_jump_expression!( + expression[area_name, t], + get_variable_multiplier(U(), d, W()), + ) + else + _add_to_jump_expression!( + expression[area_name, t], + variable[name, t], + get_variable_multiplier(U(), d, W()), + ) + end + end end return end @@ -676,7 +694,7 @@ function add_to_expression!( network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: ActivePowerBalance, - U <: OnVariable, + U <: , V <: PSY.ThermalGen, W <: Union{AbstractCompactUnitCommitment, ThermalCompactDispatch}, } From 6957939b3fb83fa2cbefc0f4e75df74882ae7ea8 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:58:46 -0600 Subject: [PATCH 17/68] remove JSON -> JSON3 --- Project.toml | 1 - src/PowerSimulations.jl | 1 - src/utils/file_utils.jl | 4 ++-- src/utils/jump_utils.jl | 2 +- 4 files changed, 3 insertions(+), 5 deletions(-) diff --git a/Project.toml b/Project.toml index 01523f9175..860ecd1727 100644 --- a/Project.toml +++ b/Project.toml @@ -38,7 +38,6 @@ DocStringExtensions = "~v0.9" HDF5 = "~0.17" InfrastructureSystems = "2" InteractiveUtils = "1" -JSON = "0.21" JSON3 = "1" JuMP = "1" LinearAlgebra = "1" diff --git a/src/PowerSimulations.jl b/src/PowerSimulations.jl index f4d882d98d..ddca095e63 100644 --- a/src/PowerSimulations.jl +++ b/src/PowerSimulations.jl @@ -412,7 +412,6 @@ import TimeSeries # I/O Imports import DataFrames -import JSON import CSV import HDF5 import PrettyTables diff --git a/src/utils/file_utils.jl b/src/utils/file_utils.jl index 9ea6116bba..8ee977d894 100644 --- a/src/utils/file_utils.jl +++ b/src/utils/file_utils.jl @@ -3,7 +3,7 @@ Return a decoded JSON file. """ function read_json(filename::AbstractString) open(filename, "r") do io - JSON.parse(io) + JSON3.parse(io) end end @@ -28,7 +28,7 @@ end function read_file_hashes(path) data = open(joinpath(path, IS.HASH_FILENAME), "r") do io - JSON.parse(io) + JSON3.parse(io) end return data["files"] diff --git a/src/utils/jump_utils.jl b/src/utils/jump_utils.jl index e610a91ae5..e8635e607e 100644 --- a/src/utils/jump_utils.jl +++ b/src/utils/jump_utils.jl @@ -6,7 +6,7 @@ function add_jump_parameter(jump_model::JuMP.Model, val::Number) end function write_data(base_power::Float64, save_path::String) - JSON.write(joinpath(save_path, "base_power.json"), JSON.json(base_power)) + JSON3.write(joinpath(save_path, "base_power.json"), JSON3.json(base_power)) return end From 8c718ce494e379826b35a9cb22c7b1b30f8097c7 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:05:55 -0600 Subject: [PATCH 18/68] change use of onvar --- src/devices_models/devices/common/add_to_expression.jl | 2 +- .../devices/common/objective_function/common.jl | 3 +++ src/devices_models/devices/common/range_constraint.jl | 4 ++-- 3 files changed, 6 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 3c7ce3277b..ca464980ae 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -694,7 +694,7 @@ function add_to_expression!( network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: ActivePowerBalance, - U <: , + U <: OnVariable, V <: PSY.ThermalGen, W <: Union{AbstractCompactUnitCommitment, ThermalCompactDispatch}, } diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index 082cb79932..3d665d35e4 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -27,6 +27,7 @@ function add_shut_down_cost!( ) where {T <: PSY.Component, U <: VariableType, V <: AbstractDeviceFormulation} multiplier = objective_function_multiplier(U(), V()) for d in devices + PSY.get_must_run(d) && continue op_cost_data = PSY.get_operation_cost(d) cost_term = shut_down_cost(op_cost_data, d, V()) iszero(cost_term) && continue @@ -114,6 +115,7 @@ function _add_start_up_cost_to_objective!( op_cost::Union{PSY.ThermalGenerationCost, PSY.MarketBidCost}, ::U, ) where {T <: VariableType, U <: AbstractDeviceFormulation} + PSY.get_must_run(component) && return cost_term = start_up_cost(op_cost, component, U()) iszero(cost_term) && return multiplier = objective_function_multiplier(T(), U()) @@ -136,6 +138,7 @@ function _add_start_up_cost_to_objective!( op_cost::PSY.ThermalGenerationCost, ::U, ) where {T <: VariableType, U <: ThermalMultiStartUnitCommitment} + PSY.get_must_run(component) && return cost_terms = start_up_cost(op_cost, component, U()) cost_term = cost_terms[MULTI_START_COST_MAP[T]] iszero(cost_term) && return diff --git a/src/devices_models/devices/common/range_constraint.jl b/src/devices_models/devices/common/range_constraint.jl index d92c5620e1..0a7c4b804c 100644 --- a/src/devices_models/devices/common/range_constraint.jl +++ b/src/devices_models/devices/common/range_constraint.jl @@ -277,7 +277,7 @@ function _add_semicontinuous_lower_bound_range_constraints_impl!( for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - if PSY.get_must_run_device(device) + if PSY.get_must_run(device) for t in time_steps con_lb[ci_name, t] = JuMP.@constraint( get_jump_model(container), @@ -311,7 +311,7 @@ function _add_semicontinuous_upper_bound_range_constraints_impl!( for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - if PSY.get_must_run_device(device) + if PSY.get_must_run(device) for t in time_steps con_ub[ci_name, t] = JuMP.@constraint( get_jump_model(container), From 0ce9578f46947dd5349ccfcedabdb1985795b764 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:03 -0600 Subject: [PATCH 19/68] add new default --- .../devices/thermal_generation.jl | 49 +++++++++++-------- 1 file changed, 29 insertions(+), 20 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 709593521e..5f61b348ef 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -61,7 +61,7 @@ get_expression_multiplier(::OnStatusParameter, ::ActivePowerRangeExpressionLB, d get_expression_multiplier(::OnStatusParameter, ::ActivePowerBalance, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power_limits(d).min #################### Initial Conditions for models ############### -initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_status(d) +initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = max(PSY.get_must_run(d), PSY.get_status(d)) initial_condition_variable(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = OnVariable() initial_condition_default(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power(d) initial_condition_variable(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = ActivePowerVariable() @@ -307,7 +307,7 @@ function add_variable!( devices::U, formulation::AbstractThermalFormulation, ) where { - T <: OnVariable, + T <: Union{OnVariable, StartVariable, StopVariable}, U <: Union{Vector{D}, IS.FlattenIteratorWrapper{D}}, } where {D <: PSY.ThermalGen} @assert !isempty(devices) @@ -719,26 +719,35 @@ function add_constraints!( for ic in initial_conditions name = PSY.get_name(get_component(ic)) - constraint[name, 1] = JuMP.@constraint( - get_jump_model(container), - varon[name, 1] == get_value(ic) + varstart[name, 1] - varstop[name, 1] - ) - aux_constraint[name, 1] = JuMP.@constraint( - get_jump_model(container), - varstart[name, 1] + varstop[name, 1] <= 1.0 - ) + if !PSY.get_must_run(get_component(ic)) + constraint[name, 1] = JuMP.@constraint( + get_jump_model(container), + varon[name, 1] == get_value(ic) + varstart[name, 1] - varstop[name, 1] + ) + aux_constraint[name, 1] = JuMP.@constraint( + get_jump_model(container), + varstart[name, 1] + varstop[name, 1] <= 1.0 + ) + end end - for t in time_steps[2:end], ic in initial_conditions - name = get_component_name(ic) - constraint[name, t] = JuMP.@constraint( - get_jump_model(container), - varon[name, t] == varon[name, t - 1] + varstart[name, t] - varstop[name, t] - ) - aux_constraint[name, t] = JuMP.@constraint( - get_jump_model(container), - varstart[name, t] + varstop[name, t] <= 1.0 - ) + for ic in initial_conditions + if PSY.get_must_run(get_component(ic)) + continue + else + name = get_component_name(ic) + for t in time_steps[2:end] + constraint[name, t] = JuMP.@constraint( + get_jump_model(container), + varon[name, t] == + varon[name, t - 1] + varstart[name, t] - varstop[name, t] + ) + aux_constraint[name, t] = JuMP.@constraint( + get_jump_model(container), + varstart[name, t] + varstop[name, t] <= 1.0 + ) + end + end end return end From 1e5011f9ede9f2970e9690a9395cf56f4c344fb5 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:17 -0600 Subject: [PATCH 20/68] WIP fix initial conditions --- src/initial_conditions/add_initial_condition.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/initial_conditions/add_initial_condition.jl b/src/initial_conditions/add_initial_condition.jl index 5dd0dba455..2d28e418d0 100644 --- a/src/initial_conditions/add_initial_condition.jl +++ b/src/initial_conditions/add_initial_condition.jl @@ -16,7 +16,7 @@ function _get_initial_conditions_value( else val = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -39,7 +39,7 @@ function _get_initial_conditions_value( else val = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -66,7 +66,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -93,7 +93,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -120,7 +120,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -147,7 +147,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -165,7 +165,7 @@ function _get_initial_conditions_value( } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -183,7 +183,7 @@ function _get_initial_conditions_value( } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end From aea5144bce72434f88b4c510eb3832e371cd0173 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:23 -0600 Subject: [PATCH 21/68] improve testing --- test/test_device_thermal_generation_constructors.jl | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index c7387981c3..c6bedc670d 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -886,7 +886,7 @@ end sys_5 = build_system(PSITestSystems, "c_sys5_uc") template_uc = ProblemTemplate(NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(sys_5))) - set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment) + set_device_model!(template_uc, ThermalStandard, ThermalStandardUnitCommitment) set_device_model!(template_uc, RenewableDispatch, FixedOutput) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) set_device_model!(template_uc, DeviceModel(Line, StaticBranch)) @@ -904,7 +904,10 @@ end solve!(model; output_dir = mktempdir()) ptdf_vars = get_variable_values(OptimizationProblemResults(model)) + power = + ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")] on = ptdf_vars[PowerSimulations.VariableKey{OnVariable, ThermalStandard}("")] - on_sundance = on[!, "Sundance"] - @test all(isapprox.(on_sundance, 1.0)) + power_sundance = power[!, "Sundance"] + @test all(isapprox.(power_sundance, 1.0)) + @test "Sundance" ∉ names(on) end From d1912be0961b435fa12ebff039f1c9c05f63caa8 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:15 -0600 Subject: [PATCH 22/68] save fixed value for variable --- src/core/optimization_container.jl | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index e3e7f2ffb7..c0fc9dfba6 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -1716,6 +1716,9 @@ function _process_duals(container::OptimizationContainer, lp_optimizer) if JuMP.has_upper_bound(first(variable)) cache[key][:ub] = JuMP.upper_bound.(variable) end + if JuMP.is_fixed(first(variable)) && is_integer_flag + cache[key][:fixed_int_value] = jump_value.(v) + end cache[key][:integer] = is_integer_flag JuMP.fix.(variable, var_cache[key]; force = true) end @@ -1756,6 +1759,9 @@ function _process_duals(container::OptimizationContainer, lp_optimizer) else JuMP.unfix.(variable) JuMP.set_binary.(variable) + if haskey(cache[key], :fixed_int_value) + JuMP.fix.(variable, cache[key][:fixed_int_value]) + end #= Needed if a model has integer variables if haskey(cache[key], :lb) && JuMP.has_lower_bound(first(variable)) JuMP.set_lower_bound.(variable, cache[key][:lb]) From 42ac5e310e0909dce023701b3c5318cd71fb0f0f Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:34 -0600 Subject: [PATCH 23/68] change implementation of must_run vars --- .../devices/thermal_generation.jl | 28 ++++++++++--------- 1 file changed, 15 insertions(+), 13 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 28b12488da..709593521e 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -319,23 +319,25 @@ function add_variable!( container, variable_type, D, - [PSY.get_name(d) for d in devices], + [PSY.get_name(d) for d in devices if !PSY.get_must_run(d)], time_steps, ) - for t in time_steps, d in devices - name = PSY.get_name(d) - variable[name, t] = JuMP.@variable( - get_jump_model(container), - base_name = "$(T)_$(D)_{$(name), $(t)}", - binary = binary - ) - if get_warm_start(settings) - init = get_variable_warm_start_value(variable_type, d, formulation) - init !== nothing && JuMP.set_start_value(variable[name, t], init) - end + for d in devices if PSY.get_must_run(d) - JuMP.fix(variable[name, t], 1.0; force = true) + continue + end + name = PSY.get_name(d) + for t in time_steps + variable[name, t] = JuMP.@variable( + get_jump_model(container), + base_name = "$(T)_$(D)_{$(name), $(t)}", + binary = binary + ) + if get_warm_start(settings) + init = get_variable_warm_start_value(variable_type, d, formulation) + init !== nothing && JuMP.set_start_value(variable[name, t], init) + end end end From ddcb5ddabc2b0879301e61b92fc3bd857b567fc3 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:50:52 -0600 Subject: [PATCH 24/68] change implementation of range constraints with must run --- .../devices/common/range_constraint.jl | 102 ++++++++++++++---- 1 file changed, 82 insertions(+), 20 deletions(-) diff --git a/src/devices_models/devices/common/range_constraint.jl b/src/devices_models/devices/common/range_constraint.jl index a4e0566066..d92c5620e1 100644 --- a/src/devices_models/devices/common/range_constraint.jl +++ b/src/devices_models/devices/common/range_constraint.jl @@ -104,7 +104,7 @@ function _add_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_lb[ci_name, t] = - JuMP.@constraint(container.JuMPmodel, array[ci_name, t] >= limits.min) + JuMP.@constraint(get_jump_model(container), array[ci_name, t] >= limits.min) end return end @@ -126,7 +126,7 @@ function _add_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_ub[ci_name, t] = - JuMP.@constraint(container.JuMPmodel, array[ci_name, t] <= limits.max) + JuMP.@constraint(get_jump_model(container), array[ci_name, t] <= limits.max) end return end @@ -246,20 +246,86 @@ function _add_semicontinuous_lower_bound_range_constraints_impl!( ) where {T <: ConstraintType, V <: PSY.Component, W <: AbstractDeviceFormulation} time_steps = get_time_steps(container) names = [PSY.get_name(d) for d in devices] - binary_variables = [OnVariable()] + con_lb = add_constraints_container!(container, T(), V, names, time_steps; meta = "lb") + varbin = get_variable(container, OnVariable(), V) + + for device in devices + ci_name = PSY.get_name(device) + limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min * varbin[ci_name, t] + ) + end + end + return +end +function _add_semicontinuous_lower_bound_range_constraints_impl!( + container::OptimizationContainer, + ::Type{T}, + array, + devices::IS.FlattenIteratorWrapper{V}, + ::DeviceModel{V, W}, +) where {T <: ConstraintType, V <: PSY.ThermalGen, W <: AbstractDeviceFormulation} + time_steps = get_time_steps(container) + names = [PSY.get_name(d) for d in devices] con_lb = add_constraints_container!(container, T(), V, names, time_steps; meta = "lb") + varbin = get_variable(container, OnVariable(), V) - @assert length(binary_variables) == 1 "Expected $(binary_variables) for $U $V $T $W to be length 1" - varbin = get_variable(container, only(binary_variables), V) + for device in devices + ci_name = PSY.get_name(device) + limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type + if PSY.get_must_run_device(device) + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min + ) + end + else + for t in time_steps + con_lb[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] >= limits.min * varbin[ci_name, t] + ) + end + end + end + return +end - for device in devices, t in time_steps +function _add_semicontinuous_upper_bound_range_constraints_impl!( + container::OptimizationContainer, + ::Type{T}, + array, + devices::IS.FlattenIteratorWrapper{V}, + model::DeviceModel{V, W}, +) where {T <: ConstraintType, V <: PSY.ThermalGen, W <: AbstractDeviceFormulation} + time_steps = get_time_steps(container) + names = [PSY.get_name(d) for d in devices] + con_ub = add_constraints_container!(container, T(), V, names, time_steps; meta = "ub") + varbin = get_variable(container, OnVariable(), V) + + for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, - array[ci_name, t] >= limits.min * varbin[ci_name, t] - ) + if PSY.get_must_run_device(device) + for t in time_steps + con_ub[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] <= limits.max + ) + end + else + for t in time_steps + con_ub[ci_name, t] = JuMP.@constraint( + get_jump_model(container), + array[ci_name, t] <= limits.max * varbin[ci_name, t] + ) + end + end end return end @@ -273,18 +339,14 @@ function _add_semicontinuous_upper_bound_range_constraints_impl!( ) where {T <: ConstraintType, V <: PSY.Component, W <: AbstractDeviceFormulation} time_steps = get_time_steps(container) names = [PSY.get_name(d) for d in devices] - binary_variables = [OnVariable()] - con_ub = add_constraints_container!(container, T(), V, names, time_steps; meta = "ub") - - @assert length(binary_variables) == 1 "Expected $(binary_variables) for $U $V $T $W to be length 1" - varbin = get_variable(container, only(binary_variables), V) + varbin = get_variable(container, OnVariable(), V) for device in devices, t in time_steps ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * varbin[ci_name, t] ) end @@ -375,7 +437,7 @@ function _add_reserve_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] >= limits.min * (1 - varbin[ci_name, t]) ) end @@ -409,7 +471,7 @@ function _add_reserve_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * (1 - varbin[ci_name, t]) ) end @@ -513,7 +575,7 @@ function _add_reserve_lower_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, X) # depends on constraint type and formulation type con_lb[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] >= limits.min * varbin[ci_name, t] ) end @@ -545,7 +607,7 @@ function _add_reserve_upper_bound_range_constraints_impl!( ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, X) # depends on constraint type and formulation type con_ub[ci_name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), array[ci_name, t] <= limits.max * varbin[ci_name, t] ) end From 0dab16735161b51e72d29e26c5e3794a53319d46 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:57:14 -0600 Subject: [PATCH 25/68] WP fix addition to expressions --- .../devices/common/add_to_expression.jl | 46 +++++++++++++------ 1 file changed, 32 insertions(+), 14 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index efbda4845c..3c7ce3277b 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -530,15 +530,24 @@ 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) - for d in devices, t in get_time_steps(container) + for d in devices name = PSY.get_name(d) bus_no_ = PSY.get_number(PSY.get_bus(d)) bus_no = PNM.get_mapped_bus_number(radial_network_reduction, bus_no_) - _add_to_jump_expression!( - expression[bus_no, t], - variable[name, t], - get_variable_multiplier(U(), d, W()), - ) + for t in get_time_steps(container) + if PSY.get_must_run(d) + _add_to_jump_expression!( + expression[bus_no, t], + get_variable_multiplier(U(), d, W()), + ) + else + _add_to_jump_expression!( + expression[bus_no, t], + variable[name, t], + get_variable_multiplier(U(), d, W()), + ) + end + end end return end @@ -558,15 +567,24 @@ function add_to_expression!( } variable = get_variable(container, U(), V) expression = get_expression(container, T(), PSY.ACBus) - for d in devices, t in get_time_steps(container) + for d in devices + name = PSY.get_name(d) bus = PSY.get_bus(d) area_name = PSY.get_name(PSY.get_area(bus)) - name = PSY.get_name(d) - _add_to_jump_expression!( - expression[area_name, t], - variable[name, t], - get_variable_multiplier(U(), d, W()), - ) + for t in get_time_steps(container) + if PSY.get_must_run(d) + _add_to_jump_expression!( + expression[area_name, t], + get_variable_multiplier(U(), d, W()), + ) + else + _add_to_jump_expression!( + expression[area_name, t], + variable[name, t], + get_variable_multiplier(U(), d, W()), + ) + end + end end return end @@ -676,7 +694,7 @@ function add_to_expression!( network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: ActivePowerBalance, - U <: OnVariable, + U <: , V <: PSY.ThermalGen, W <: Union{AbstractCompactUnitCommitment, ThermalCompactDispatch}, } From 36d176508c8e30cc9b01d1dede43844719e817b7 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 20 Sep 2024 20:58:46 -0600 Subject: [PATCH 26/68] remove JSON -> JSON3 --- Project.toml | 1 - src/PowerSimulations.jl | 1 - src/utils/file_utils.jl | 4 ++-- src/utils/jump_utils.jl | 2 +- 4 files changed, 3 insertions(+), 5 deletions(-) diff --git a/Project.toml b/Project.toml index 53a901c19e..e2c8fafa10 100644 --- a/Project.toml +++ b/Project.toml @@ -38,7 +38,6 @@ DocStringExtensions = "~v0.9" HDF5 = "~0.17" InfrastructureSystems = "2" InteractiveUtils = "1" -JSON = "0.21" JSON3 = "1" JuMP = "1" LinearAlgebra = "1" diff --git a/src/PowerSimulations.jl b/src/PowerSimulations.jl index f4d882d98d..ddca095e63 100644 --- a/src/PowerSimulations.jl +++ b/src/PowerSimulations.jl @@ -412,7 +412,6 @@ import TimeSeries # I/O Imports import DataFrames -import JSON import CSV import HDF5 import PrettyTables diff --git a/src/utils/file_utils.jl b/src/utils/file_utils.jl index 9ea6116bba..8ee977d894 100644 --- a/src/utils/file_utils.jl +++ b/src/utils/file_utils.jl @@ -3,7 +3,7 @@ Return a decoded JSON file. """ function read_json(filename::AbstractString) open(filename, "r") do io - JSON.parse(io) + JSON3.parse(io) end end @@ -28,7 +28,7 @@ end function read_file_hashes(path) data = open(joinpath(path, IS.HASH_FILENAME), "r") do io - JSON.parse(io) + JSON3.parse(io) end return data["files"] diff --git a/src/utils/jump_utils.jl b/src/utils/jump_utils.jl index e610a91ae5..e8635e607e 100644 --- a/src/utils/jump_utils.jl +++ b/src/utils/jump_utils.jl @@ -6,7 +6,7 @@ function add_jump_parameter(jump_model::JuMP.Model, val::Number) end function write_data(base_power::Float64, save_path::String) - JSON.write(joinpath(save_path, "base_power.json"), JSON.json(base_power)) + JSON3.write(joinpath(save_path, "base_power.json"), JSON3.json(base_power)) return end From 236eda625ab9e4927beebfe0473ba616fbb7fff2 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:05:55 -0600 Subject: [PATCH 27/68] change use of onvar --- src/devices_models/devices/common/add_to_expression.jl | 2 +- .../devices/common/objective_function/common.jl | 3 +++ src/devices_models/devices/common/range_constraint.jl | 4 ++-- 3 files changed, 6 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 3c7ce3277b..ca464980ae 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -694,7 +694,7 @@ function add_to_expression!( network_model::NetworkModel{CopperPlatePowerModel}, ) where { T <: ActivePowerBalance, - U <: , + U <: OnVariable, V <: PSY.ThermalGen, W <: Union{AbstractCompactUnitCommitment, ThermalCompactDispatch}, } diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index 082cb79932..3d665d35e4 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -27,6 +27,7 @@ function add_shut_down_cost!( ) where {T <: PSY.Component, U <: VariableType, V <: AbstractDeviceFormulation} multiplier = objective_function_multiplier(U(), V()) for d in devices + PSY.get_must_run(d) && continue op_cost_data = PSY.get_operation_cost(d) cost_term = shut_down_cost(op_cost_data, d, V()) iszero(cost_term) && continue @@ -114,6 +115,7 @@ function _add_start_up_cost_to_objective!( op_cost::Union{PSY.ThermalGenerationCost, PSY.MarketBidCost}, ::U, ) where {T <: VariableType, U <: AbstractDeviceFormulation} + PSY.get_must_run(component) && return cost_term = start_up_cost(op_cost, component, U()) iszero(cost_term) && return multiplier = objective_function_multiplier(T(), U()) @@ -136,6 +138,7 @@ function _add_start_up_cost_to_objective!( op_cost::PSY.ThermalGenerationCost, ::U, ) where {T <: VariableType, U <: ThermalMultiStartUnitCommitment} + PSY.get_must_run(component) && return cost_terms = start_up_cost(op_cost, component, U()) cost_term = cost_terms[MULTI_START_COST_MAP[T]] iszero(cost_term) && return diff --git a/src/devices_models/devices/common/range_constraint.jl b/src/devices_models/devices/common/range_constraint.jl index d92c5620e1..0a7c4b804c 100644 --- a/src/devices_models/devices/common/range_constraint.jl +++ b/src/devices_models/devices/common/range_constraint.jl @@ -277,7 +277,7 @@ function _add_semicontinuous_lower_bound_range_constraints_impl!( for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - if PSY.get_must_run_device(device) + if PSY.get_must_run(device) for t in time_steps con_lb[ci_name, t] = JuMP.@constraint( get_jump_model(container), @@ -311,7 +311,7 @@ function _add_semicontinuous_upper_bound_range_constraints_impl!( for device in devices ci_name = PSY.get_name(device) limits = get_min_max_limits(device, T, W) # depends on constraint type and formulation type - if PSY.get_must_run_device(device) + if PSY.get_must_run(device) for t in time_steps con_ub[ci_name, t] = JuMP.@constraint( get_jump_model(container), From c734e780fb7f8930bd0ecda908a0a62a64d7453e Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:03 -0600 Subject: [PATCH 28/68] add new default --- .../devices/thermal_generation.jl | 49 +++++++++++-------- 1 file changed, 29 insertions(+), 20 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 709593521e..5f61b348ef 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -61,7 +61,7 @@ get_expression_multiplier(::OnStatusParameter, ::ActivePowerRangeExpressionLB, d get_expression_multiplier(::OnStatusParameter, ::ActivePowerBalance, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power_limits(d).min #################### Initial Conditions for models ############### -initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_status(d) +initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = max(PSY.get_must_run(d), PSY.get_status(d)) initial_condition_variable(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = OnVariable() initial_condition_default(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power(d) initial_condition_variable(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = ActivePowerVariable() @@ -307,7 +307,7 @@ function add_variable!( devices::U, formulation::AbstractThermalFormulation, ) where { - T <: OnVariable, + T <: Union{OnVariable, StartVariable, StopVariable}, U <: Union{Vector{D}, IS.FlattenIteratorWrapper{D}}, } where {D <: PSY.ThermalGen} @assert !isempty(devices) @@ -719,26 +719,35 @@ function add_constraints!( for ic in initial_conditions name = PSY.get_name(get_component(ic)) - constraint[name, 1] = JuMP.@constraint( - get_jump_model(container), - varon[name, 1] == get_value(ic) + varstart[name, 1] - varstop[name, 1] - ) - aux_constraint[name, 1] = JuMP.@constraint( - get_jump_model(container), - varstart[name, 1] + varstop[name, 1] <= 1.0 - ) + if !PSY.get_must_run(get_component(ic)) + constraint[name, 1] = JuMP.@constraint( + get_jump_model(container), + varon[name, 1] == get_value(ic) + varstart[name, 1] - varstop[name, 1] + ) + aux_constraint[name, 1] = JuMP.@constraint( + get_jump_model(container), + varstart[name, 1] + varstop[name, 1] <= 1.0 + ) + end end - for t in time_steps[2:end], ic in initial_conditions - name = get_component_name(ic) - constraint[name, t] = JuMP.@constraint( - get_jump_model(container), - varon[name, t] == varon[name, t - 1] + varstart[name, t] - varstop[name, t] - ) - aux_constraint[name, t] = JuMP.@constraint( - get_jump_model(container), - varstart[name, t] + varstop[name, t] <= 1.0 - ) + for ic in initial_conditions + if PSY.get_must_run(get_component(ic)) + continue + else + name = get_component_name(ic) + for t in time_steps[2:end] + constraint[name, t] = JuMP.@constraint( + get_jump_model(container), + varon[name, t] == + varon[name, t - 1] + varstart[name, t] - varstop[name, t] + ) + aux_constraint[name, t] = JuMP.@constraint( + get_jump_model(container), + varstart[name, t] + varstop[name, t] <= 1.0 + ) + end + end end return end From 8ce9292e3db638cc8d433a48c8cc41d256f4fd25 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:17 -0600 Subject: [PATCH 29/68] WIP fix initial conditions --- src/initial_conditions/add_initial_condition.jl | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/initial_conditions/add_initial_condition.jl b/src/initial_conditions/add_initial_condition.jl index 5dd0dba455..2d28e418d0 100644 --- a/src/initial_conditions/add_initial_condition.jl +++ b/src/initial_conditions/add_initial_condition.jl @@ -16,7 +16,7 @@ function _get_initial_conditions_value( else val = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -39,7 +39,7 @@ function _get_initial_conditions_value( else val = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -66,7 +66,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -93,7 +93,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -120,7 +120,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end @@ -147,7 +147,7 @@ function _get_initial_conditions_value( val = PSY.get_time_at_status(component) end end - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -165,7 +165,7 @@ function _get_initial_conditions_value( } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, add_jump_parameter(get_jump_model(container), val)) end @@ -183,7 +183,7 @@ function _get_initial_conditions_value( } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) - @debug "Device $(PSY.get_name(component)) initialized DeviceStatus as $var_type" _group = + @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS return T(component, val) end From 0b6baac44c907058c12144c0307f5fed9b556a44 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 13:06:23 -0600 Subject: [PATCH 30/68] improve testing --- test/test_device_thermal_generation_constructors.jl | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index c7387981c3..c6bedc670d 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -886,7 +886,7 @@ end sys_5 = build_system(PSITestSystems, "c_sys5_uc") template_uc = ProblemTemplate(NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(sys_5))) - set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment) + set_device_model!(template_uc, ThermalStandard, ThermalStandardUnitCommitment) set_device_model!(template_uc, RenewableDispatch, FixedOutput) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) set_device_model!(template_uc, DeviceModel(Line, StaticBranch)) @@ -904,7 +904,10 @@ end solve!(model; output_dir = mktempdir()) ptdf_vars = get_variable_values(OptimizationProblemResults(model)) + power = + ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")] on = ptdf_vars[PowerSimulations.VariableKey{OnVariable, ThermalStandard}("")] - on_sundance = on[!, "Sundance"] - @test all(isapprox.(on_sundance, 1.0)) + power_sundance = power[!, "Sundance"] + @test all(isapprox.(power_sundance, 1.0)) + @test "Sundance" ∉ names(on) end From 2204b5bc8515d6f7601c1102be25322b8e68c7f4 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 23 Sep 2024 14:20:18 -0600 Subject: [PATCH 31/68] add lp writter option --- src/operation/operation_model_interface.jl | 1 + src/utils/jump_utils.jl | 7 +++++++ 2 files changed, 8 insertions(+) diff --git a/src/operation/operation_model_interface.jl b/src/operation/operation_model_interface.jl index d2157ee967..704331b09a 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -114,6 +114,7 @@ function solve_impl!(model::OperationModel) tss = replace("$(ts)", ":" => "_") model_export_path = joinpath(model_output_dir, "exported_$(model_name)_$(tss).json") serialize_optimization_model(container, model_export_path) + @debug write_lp_file(get_jump_model(container), replace(model_export_path, ".json" => ".lp")) end status = solve_impl!(container, get_system(model)) diff --git a/src/utils/jump_utils.jl b/src/utils/jump_utils.jl index e8635e607e..129414c0d9 100644 --- a/src/utils/jump_utils.jl +++ b/src/utils/jump_utils.jl @@ -315,6 +315,13 @@ function serialize_jump_optimization_model(jump_model::JuMP.Model, save_path::St return end +function write_lp_file(jump_model::JuMP.Model, save_path::String) + MOF_model = MOPFM(; format = MOI.FileFormats.FORMAT_LP) + MOI.copy_to(MOF_model, JuMP.backend(jump_model)) + MOI.write_to_file(MOF_model, save_path) + return +end + # check_conflict_status functions can't be tested on CI because free solvers don't support IIS function check_conflict_status( jump_model::JuMP.Model, From 066ac2d7f66c8f0f672d837187032a95ba42ddbd Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 24 Sep 2024 13:23:00 -0600 Subject: [PATCH 32/68] WIP change IC handling for must run --- src/core/initial_conditions.jl | 2 +- src/core/optimization_container.jl | 3 +- .../device_constructors/branch_constructor.jl | 2 +- .../devices/common/duration_constraints.jl | 41 ++++++--- .../devices/thermal_generation.jl | 16 ++-- .../add_initial_condition.jl | 91 ++++++++++++++----- src/operation/operation_model_interface.jl | 5 +- ..._device_thermal_generation_constructors.jl | 11 ++- 8 files changed, 118 insertions(+), 53 deletions(-) diff --git a/src/core/initial_conditions.jl b/src/core/initial_conditions.jl index 329cb47c94..5474275422 100644 --- a/src/core/initial_conditions.jl +++ b/src/core/initial_conditions.jl @@ -3,7 +3,7 @@ Container for the initial condition data """ mutable struct InitialCondition{ T <: InitialConditionType, - U <: Union{JuMP.VariableRef, Float64}, + U <: Union{JuMP.VariableRef, Float64, Nothing}, } component::PSY.Component value::U diff --git a/src/core/optimization_container.jl b/src/core/optimization_container.jl index c0fc9dfba6..084bd32103 100644 --- a/src/core/optimization_container.jl +++ b/src/core/optimization_container.jl @@ -1503,7 +1503,8 @@ function _add_initial_condition_container!( else param_type = Float64 end - ini_conds = Vector{InitialCondition{T, param_type}}(undef, length_devices) + ini_type = Union{InitialCondition{T, param_type}, InitialCondition{T, Nothing}} + ini_conds = Vector{ini_type}(undef, length_devices) _assign_container!(container.initial_conditions, ic_key, ini_conds) return ini_conds end diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 9aff973682..616e82d74b 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -167,7 +167,7 @@ function construct_device!( end function construct_device!( - ::OptimizationContainer, + container::OptimizationContainer, sys::PSY.System, ::ModelConstructStage, model::DeviceModel{<:PSY.ACBranch, StaticBranchUnbounded}, diff --git a/src/devices_models/devices/common/duration_constraints.jl b/src/devices_models/devices/common/duration_constraints.jl index 954a7e2aa5..12a1f8d7b8 100644 --- a/src/devices_models/devices/common/duration_constraints.jl +++ b/src/devices_models/devices/common/duration_constraints.jl @@ -52,7 +52,10 @@ function device_duration_retrospective!( varstart = get_variable(container, var_types[2], T) varstop = get_variable(container, var_types[3], T) - set_names = [get_component_name(ic) for ic in initial_duration[:, 1]] + set_names = [ + get_component_name(ic) for + ic in initial_duration[:, 1] if !isnothing(get_value(ic)) + ] con_up = add_constraints_container!( container, cons_type, @@ -72,6 +75,7 @@ function device_duration_retrospective!( for t in time_steps for (ix, ic) in enumerate(initial_duration[:, 1]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Up-time Constraint lhs_on = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -84,10 +88,11 @@ function device_duration_retrospective!( JuMP.add_to_expression!(lhs_on, 1) end con_up[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_on - varon[name, t] <= 0.0) + JuMP.@constraint(get_jump_model(container), lhs_on - varon[name, t] <= 0.0) end for (ix, ic) in enumerate(initial_duration[:, 2]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Down-time Constraint lhs_off = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -100,11 +105,12 @@ function device_duration_retrospective!( JuMP.add_to_expression!(lhs_off, 1) end con_down[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_off + varon[name, t] <= 1.0) + JuMP.@constraint(get_jump_model(container), lhs_off + varon[name, t] <= 1.0) end end return end + @doc raw""" This formulation of the duration constraints looks ahead in the time frame of the model. @@ -165,6 +171,7 @@ function device_duration_look_ahead!( for t in time_steps for (ix, ic) in enumerate(initial_duration[:, 1]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Up-time Constraint lhs_on = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -177,12 +184,13 @@ function device_duration_look_ahead!( lhs_on += get_value(ic) end con_up[name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), varstop[name, t] * duration_data[ix].up - lhs_on <= 0.0 ) end for (ix, ic) in enumerate(initial_duration[:, 2]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Down-time Constraint lhs_off = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -195,7 +203,7 @@ function device_duration_look_ahead!( lhs_off += get_value(ic) end con_down[name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), varstart[name, t] * duration_data[ix].down - lhs_off <= 0.0 ) end @@ -278,8 +286,8 @@ function device_duration_parameters!( for t in time_steps for (ix, ic) in enumerate(initial_duration[:, 1]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) - # Minimum Up-time Constraint lhs_on = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) for i in UnitRange{Int}(Int(t - duration_data[ix].up + 1), t) @@ -294,16 +302,20 @@ function device_duration_parameters!( if t <= duration_data[ix].up lhs_on += get_value(ic) con_up[name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), varstop[name, t] * duration_data[ix].up - lhs_on <= 0.0 ) else con_up[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_on - varon[name, t] <= 0.0) + JuMP.@constraint( + get_jump_model(container), + lhs_on - varon[name, t] <= 0.0 + ) end end for (ix, ic) in enumerate(initial_duration[:, 2]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Down-time Constraint lhs_off = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -319,12 +331,15 @@ function device_duration_parameters!( if t <= duration_data[ix].down lhs_off += get_value(ic) con_down[name, t] = JuMP.@constraint( - container.JuMPmodel, + get_jump_model(container), varstart[name, t] * duration_data[ix].down - lhs_off <= 0.0 ) else con_down[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_off + varon[name, t] <= 1.0) + JuMP.@constraint( + get_jump_model(container), + lhs_off + varon[name, t] <= 1.0 + ) end end end @@ -395,6 +410,7 @@ function device_duration_compact_retrospective!( total_time_steps = length(time_steps) for t in time_steps for (ix, ic) in enumerate(initial_duration[:, 1]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Up-time Constraint lhs_on = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -413,10 +429,11 @@ function device_duration_compact_retrospective!( continue end con_up[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_on - varon[name, t] <= 0.0) + JuMP.@constraint(get_jump_model(container), lhs_on - varon[name, t] <= 0.0) end for (ix, ic) in enumerate(initial_duration[:, 2]) + isnothing(get_value(ic)) && continue name = get_component_name(ic) # Minimum Down-time Constraint lhs_off = JuMP.GenericAffExpr{Float64, JuMP.VariableRef}(0) @@ -435,7 +452,7 @@ function device_duration_compact_retrospective!( continue end con_down[name, t] = - JuMP.@constraint(container.JuMPmodel, lhs_off + varon[name, t] <= 1.0) + JuMP.@constraint(get_jump_model(container), lhs_off + varon[name, t] <= 1.0) end end for c in [con_up, con_down] diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 5f61b348ef..bf8f5cbc7c 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -61,7 +61,7 @@ get_expression_multiplier(::OnStatusParameter, ::ActivePowerRangeExpressionLB, d get_expression_multiplier(::OnStatusParameter, ::ActivePowerBalance, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power_limits(d).min #################### Initial Conditions for models ############### -initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = max(PSY.get_must_run(d), PSY.get_status(d)) +initial_condition_default(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_status(d) ? 1.0 : 0.0 initial_condition_variable(::DeviceStatus, d::PSY.ThermalGen, ::AbstractThermalFormulation) = OnVariable() initial_condition_default(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_active_power(d) initial_condition_variable(::DevicePower, d::PSY.ThermalGen, ::AbstractThermalFormulation) = ActivePowerVariable() @@ -818,14 +818,14 @@ function calculate_aux_variable_value!( time_steps = get_time_steps(container) for ix in eachindex(JuMP.axes(aux_variable_container)[1]) - IS.@assert_op JuMP.axes(aux_variable_container)[1][ix] == - JuMP.axes(on_variable_results)[1][ix] - IS.@assert_op JuMP.axes(aux_variable_container)[1][ix] == - get_component_name(ini_cond[ix]) - on_var = jump_value.(on_variable_results.data[ix, :]) ini_cond_value = get_condition(ini_cond[ix]) - aux_variable_container.data[ix, :] .= ini_cond_value - sum_on_var = sum(on_var) + if isnothing(get_value(ini_cond[ix])) + sum_on_var = time_steps[end] + else + on_var = jump_value.(on_variable_results.data[ix, :]) + aux_variable_container.data[ix, :] .= ini_cond_value + sum_on_var = sum(on_var) + end if sum_on_var == time_steps[end] # Unit was always on aux_variable_container.data[ix, :] += time_steps elseif sum_on_var == 0.0 # Unit was always off diff --git a/src/initial_conditions/add_initial_condition.jl b/src/initial_conditions/add_initial_condition.jl index 2d28e418d0..3a079c572c 100644 --- a/src/initial_conditions/add_initial_condition.jl +++ b/src/initial_conditions/add_initial_condition.jl @@ -5,7 +5,7 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, Float64}, + T <: Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}, V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, W <: PSY.Component, } where {U <: InitialConditionType} @@ -18,7 +18,7 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, val) + return InitialCondition{U, Float64}(component, val) end function _get_initial_conditions_value( @@ -28,8 +28,8 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, JuMP.VariableRef}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, + T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, + V <: AbstractThermalFormulation, W <: PSY.Component, } where {U <: InitialConditionType} ic_data = get_initial_conditions_data(container) @@ -41,7 +41,10 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, add_jump_parameter(get_jump_model(container), val)) + return InitialCondition{U, JuMP.VariableRef}( + component, + add_jump_parameter(get_jump_model(container), val), + ) end function _get_initial_conditions_value( @@ -51,14 +54,14 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, Float64}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, + T <: Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}, + V <: AbstractThermalFormulation, W <: PSY.Component, } where {U <: InitialTimeDurationOff} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) if !has_initial_condition_value(ic_data, var_type, W) - val = initial_condition_default(U(), component, V()) + @show val = initial_condition_default(U(), component, V()) else var = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] val = 0.0 @@ -68,7 +71,7 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, val) + return InitialCondition{U, Float64}(component, val) end function _get_initial_conditions_value( @@ -78,9 +81,9 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, JuMP.VariableRef}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, - W <: PSY.Component, + T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, + V <: AbstractThermalFormulation, + W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOff} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) @@ -95,7 +98,10 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, add_jump_parameter(get_jump_model(container), val)) + return InitialCondition{U, JuMP.VariableRef}( + component, + add_jump_parameter(get_jump_model(container), val), + ) end function _get_initial_conditions_value( @@ -105,14 +111,14 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, Float64}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, - W <: PSY.Component, + T <: Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}, + V <: AbstractDeviceFormulation, + W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOn} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) if !has_initial_condition_value(ic_data, var_type, W) - val = initial_condition_default(U(), component, V()) + @show val = initial_condition_default(U(), component, V()) else var = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] val = 0.0 @@ -122,7 +128,7 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, val) + return InitialCondition{U, Float64}(component, val) end function _get_initial_conditions_value( @@ -132,9 +138,9 @@ function _get_initial_conditions_value( ::V, container::OptimizationContainer, ) where { - T <: InitialCondition{U, JuMP.VariableRef}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, - W <: PSY.Component, + T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, + V <: AbstractDeviceFormulation, + W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOn} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) @@ -149,7 +155,10 @@ function _get_initial_conditions_value( end @debug "Device $(PSY.get_name(component)) initialized $U as $var_type" _group = LOG_GROUP_BUILD_INITIAL_CONDITIONS - return T(component, add_jump_parameter(get_jump_model(container), val)) + return InitialCondition{U, JuMP.VariableRef}( + component, + add_jump_parameter(get_jump_model(container), val), + ) end function _get_initial_conditions_value( @@ -160,8 +169,8 @@ function _get_initial_conditions_value( container::OptimizationContainer, ) where { T <: InitialCondition{U, JuMP.VariableRef}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, - W <: PSY.Component, + V <: AbstractDeviceFormulation, + W <: PSY.ThermalGen, } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) @@ -178,7 +187,7 @@ function _get_initial_conditions_value( container::OptimizationContainer, ) where { T <: InitialCondition{U, Float64}, - V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, + V <: AbstractDeviceFormulation, W <: PSY.Component, } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) @@ -209,3 +218,35 @@ function add_initial_condition!( end return end + +function add_initial_condition!( + container::OptimizationContainer, + components::Union{Vector{T}, IS.FlattenIteratorWrapper{T}}, + ::U, + ::D, +) where { + T <: PSY.ThermalGen, + U <: AbstractThermalFormulation, + D <: InitialConditionType, +} + if get_rebuild_model(get_settings(container)) && has_container_key(container, D, T) + return + end + + ini_cond_vector = add_initial_condition_container!(container, D(), T, components) + for (ix, component) in enumerate(components) + if PSY.get_must_run(component) + ini_cond_vector[ix] = InitialCondition{D, Nothing}(component, nothing) + else + ini_cond_vector[ix] = + _get_initial_conditions_value( + ini_cond_vector, + component, + D(), + U(), + container, + ) + end + end + return +end diff --git a/src/operation/operation_model_interface.jl b/src/operation/operation_model_interface.jl index 704331b09a..e092cea01b 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -114,7 +114,10 @@ function solve_impl!(model::OperationModel) tss = replace("$(ts)", ":" => "_") model_export_path = joinpath(model_output_dir, "exported_$(model_name)_$(tss).json") serialize_optimization_model(container, model_export_path) - @debug write_lp_file(get_jump_model(container), replace(model_export_path, ".json" => ".lp")) + @debug write_lp_file( + get_jump_model(container), + replace(model_export_path, ".json" => ".lp"), + ) end status = solve_impl!(container, get_system(model)) diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index c6bedc670d..2e092a239c 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -885,24 +885,27 @@ end @testset "Test Must Run ThermalGen" begin sys_5 = build_system(PSITestSystems, "c_sys5_uc") template_uc = - ProblemTemplate(NetworkModel(PTDFPowerModel; PTDF_matrix = PTDF(sys_5))) + ProblemTemplate(NetworkModel(CopperPlatePowerModel)) set_device_model!(template_uc, ThermalStandard, ThermalStandardUnitCommitment) - set_device_model!(template_uc, RenewableDispatch, FixedOutput) + #set_device_model!(template_uc, RenewableDispatch, FixedOutput) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) - set_device_model!(template_uc, DeviceModel(Line, StaticBranch)) + set_device_model!(template_uc, DeviceModel(Line, StaticBranchUnbounded)) # Set Must Run the most expensive one: Sundance sundance = get_component(ThermalStandard, sys_5, "Sundance") + set_status!(sundance, true) set_must_run!(sundance, true) model = DecisionModel( template_uc, sys_5; name = "UC", - optimizer = HiGHS_optimizer, + optimizer = Xpress.Optimizer, system_to_file = false, + store_variable_names = true, ) solve!(model; output_dir = mktempdir()) + serialize_optimization_model(model, "branch.json") ptdf_vars = get_variable_values(OptimizationProblemResults(model)) power = ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")] From fa4c2fcd20a680a4ae3a7ce8b5d41ad2673454f1 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 24 Sep 2024 18:15:35 -0600 Subject: [PATCH 33/68] update tests --- ..._device_thermal_generation_constructors.jl | 44 +++++++++++-------- 1 file changed, 26 insertions(+), 18 deletions(-) diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index 2e092a239c..1aa73ee2af 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -893,24 +893,32 @@ end # Set Must Run the most expensive one: Sundance sundance = get_component(ThermalStandard, sys_5, "Sundance") - set_status!(sundance, true) set_must_run!(sundance, true) - model = DecisionModel( - template_uc, - sys_5; - name = "UC", - optimizer = Xpress.Optimizer, - system_to_file = false, - store_variable_names = true, - ) + for rebuild in [true, false] + model = DecisionModel( + template_uc, + sys_5; + name = "UC", + optimizer = HiGHS_optimizer, + system_to_file = false, + store_variable_names = true, + rebuild_model = rebuild, + ) - solve!(model; output_dir = mktempdir()) - serialize_optimization_model(model, "branch.json") - ptdf_vars = get_variable_values(OptimizationProblemResults(model)) - power = - ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")] - on = ptdf_vars[PowerSimulations.VariableKey{OnVariable, ThermalStandard}("")] - power_sundance = power[!, "Sundance"] - @test all(isapprox.(power_sundance, 1.0)) - @test "Sundance" ∉ names(on) + solve!(model; output_dir = mktempdir()) + serialize_optimization_model(model, "branch.json") + ptdf_vars = get_variable_values(OptimizationProblemResults(model)) + power = + ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}( + "", + )] + on = ptdf_vars[PowerSimulations.VariableKey{OnVariable, ThermalStandard}("")] + start = ptdf_vars[PowerSimulations.VariableKey{StartVariable, ThermalStandard}("")] + stop = ptdf_vars[PowerSimulations.VariableKey{StopVariable, ThermalStandard}("")] + power_sundance = power[!, "Sundance"] + @test all(power_sundance .>= 1.0) + for v in [on, start, stop] + @test "Sundance" ∉ names(v) + end + end end From a4fea254ef94893c1d34bc708b380c650a7309bf Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 24 Sep 2024 18:15:46 -0600 Subject: [PATCH 34/68] update aux vars calculations --- src/devices_models/devices/thermal_generation.jl | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index bf8f5cbc7c..bfd514ca80 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -818,10 +818,13 @@ function calculate_aux_variable_value!( time_steps = get_time_steps(container) for ix in eachindex(JuMP.axes(aux_variable_container)[1]) - ini_cond_value = get_condition(ini_cond[ix]) + # if its nothing it means the thermal unit was on must run + # so there is nothing to do but to add the total number of time steps + # to the count if isnothing(get_value(ini_cond[ix])) sum_on_var = time_steps[end] else + ini_cond_value = get_condition(ini_cond[ix]) on_var = jump_value.(on_variable_results.data[ix, :]) aux_variable_container.data[ix, :] .= ini_cond_value sum_on_var = sum(on_var) @@ -859,10 +862,11 @@ function calculate_aux_variable_value!( time_steps = get_time_steps(container) for ix in eachindex(JuMP.axes(aux_variable_container)[1]) - IS.@assert_op JuMP.axes(aux_variable_container)[1][ix] == - JuMP.axes(on_variable_results)[1][ix] - IS.@assert_op JuMP.axes(aux_variable_container)[1][ix] == - get_component_name(ini_cond[ix]) + # if its nothing it means the thermal unit was on must run + # so there is nothing to do but continue + if isnothing(get_value(ini_cond[ix])) + continue + end on_var = jump_value.(on_variable_results.data[ix, :]) ini_cond_value = get_condition(ini_cond[ix]) aux_variable_container.data[ix, :] .= ini_cond_value From 6db1560db4ca120d1be25d1382fa8763f748193d Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 24 Sep 2024 18:15:57 -0600 Subject: [PATCH 35/68] add missing method to avoid ambiguity --- .../add_initial_condition.jl | 38 ++++++++++++------- 1 file changed, 24 insertions(+), 14 deletions(-) diff --git a/src/initial_conditions/add_initial_condition.jl b/src/initial_conditions/add_initial_condition.jl index 3a079c572c..6c171a7763 100644 --- a/src/initial_conditions/add_initial_condition.jl +++ b/src/initial_conditions/add_initial_condition.jl @@ -1,3 +1,17 @@ +function _get_initial_conditions_value( + ::Vector{T}, + component::W, + ::U, + ::V, + container::OptimizationContainer, +) where { + T <: InitialCondition{U, Nothing}, + V <: Union{AbstractDeviceFormulation, AbstractServiceFormulation}, + W <: PSY.Component, +} where {U <: InitialConditionType} + return InitialCondition{U, Nothing}(component, nothing) +end + function _get_initial_conditions_value( ::Vector{T}, component::W, @@ -29,7 +43,7 @@ function _get_initial_conditions_value( container::OptimizationContainer, ) where { T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, - V <: AbstractThermalFormulation, + V <: AbstractDeviceFormulation, W <: PSY.Component, } where {U <: InitialConditionType} ic_data = get_initial_conditions_data(container) @@ -48,20 +62,19 @@ function _get_initial_conditions_value( end function _get_initial_conditions_value( - ::Vector{T}, + ::Vector{Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}}, component::W, ::U, ::V, container::OptimizationContainer, ) where { - T <: Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}, V <: AbstractThermalFormulation, W <: PSY.Component, } where {U <: InitialTimeDurationOff} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) if !has_initial_condition_value(ic_data, var_type, W) - @show val = initial_condition_default(U(), component, V()) + val = initial_condition_default(U(), component, V()) else var = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] val = 0.0 @@ -75,13 +88,12 @@ function _get_initial_conditions_value( end function _get_initial_conditions_value( - ::Vector{T}, + ::Vector{Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}}, component::W, ::U, ::V, container::OptimizationContainer, ) where { - T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, V <: AbstractThermalFormulation, W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOff} @@ -105,20 +117,19 @@ function _get_initial_conditions_value( end function _get_initial_conditions_value( - ::Vector{T}, + ::Vector{Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}}, component::W, ::U, ::V, container::OptimizationContainer, ) where { - T <: Union{InitialCondition{U, Float64}, InitialCondition{U, Nothing}}, - V <: AbstractDeviceFormulation, + V <: AbstractThermalFormulation, W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOn} ic_data = get_initial_conditions_data(container) var_type = initial_condition_variable(U(), component, V()) if !has_initial_condition_value(ic_data, var_type, W) - @show val = initial_condition_default(U(), component, V()) + val = initial_condition_default(U(), component, V()) else var = get_initial_condition_value(ic_data, var_type, W)[1, PSY.get_name(component)] val = 0.0 @@ -132,14 +143,13 @@ function _get_initial_conditions_value( end function _get_initial_conditions_value( - ::Vector{T}, + ::Vector{Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}}, component::W, ::U, ::V, container::OptimizationContainer, ) where { - T <: Union{InitialCondition{U, JuMP.VariableRef}, InitialCondition{U, Nothing}}, - V <: AbstractDeviceFormulation, + V <: AbstractThermalFormulation, W <: PSY.ThermalGen, } where {U <: InitialTimeDurationOn} ic_data = get_initial_conditions_data(container) @@ -170,7 +180,7 @@ function _get_initial_conditions_value( ) where { T <: InitialCondition{U, JuMP.VariableRef}, V <: AbstractDeviceFormulation, - W <: PSY.ThermalGen, + W <: PSY.Component, } where {U <: InitialEnergyLevel} var_type = initial_condition_variable(U(), component, V()) val = initial_condition_default(U(), component, V()) From 283b6d0cc52d49b0f94d976b6b8136e61aaed99c Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Tue, 24 Sep 2024 18:16:17 -0600 Subject: [PATCH 36/68] use correct time in inicond --- ...itial_conditions_update_in_memory_store.jl | 68 ++++++++++++------- 1 file changed, 44 insertions(+), 24 deletions(-) diff --git a/src/operation/initial_conditions_update_in_memory_store.jl b/src/operation/initial_conditions_update_in_memory_store.jl index 7c7793050b..48d56e9afd 100644 --- a/src/operation/initial_conditions_update_in_memory_store.jl +++ b/src/operation/initial_conditions_update_in_memory_store.jl @@ -2,12 +2,15 @@ ################## ic updates from store for emulation problems simulation ################# function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{InitialTimeDurationOn, Nothing}, + InitialCondition{InitialTimeDurationOn, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{InitialTimeDurationOn, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, TimeDurationOn(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -16,12 +19,15 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{InitialTimeDurationOff, Nothing}, + InitialCondition{InitialTimeDurationOff, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{InitialTimeDurationOff, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, TimeDurationOff(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -30,12 +36,15 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{DevicePower, Nothing}, + InitialCondition{DevicePower, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{DevicePower, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, ActivePowerVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -44,12 +53,15 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{DeviceStatus, Nothing}, + InitialCondition{DeviceStatus, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{DeviceStatus, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, OnVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -58,12 +70,15 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{DeviceAboveMinPower, Nothing}, + InitialCondition{DeviceAboveMinPower, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{DeviceAboveMinPower, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, PowerAboveMinimumVariable(), get_component_type(ic)) @@ -72,6 +87,7 @@ function update_initial_conditions!( return end +#= Unused without the AGC model enabled function update_initial_conditions!( ics::Vector{T}, store::EmulationModelStore, @@ -85,14 +101,18 @@ function update_initial_conditions!( end return end +=# function update_initial_conditions!( - ics::Vector{T}, + ics::Vector{ + Union{ + InitialCondition{InitialEnergyLevel, Nothing}, + InitialCondition{InitialEnergyLevel, JuMP.VariableRef}, + }, + }, store::EmulationModelStore, ::Dates.Millisecond, -) where { - T <: InitialCondition{InitialEnergyLevel, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} +) for ic in ics var_val = get_value(store, EnergyVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) From 3e2354c5172307e270677e380bd9b5254a671d35 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 11:57:48 -0600 Subject: [PATCH 37/68] remove file write --- test/test_device_thermal_generation_constructors.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/test/test_device_thermal_generation_constructors.jl b/test/test_device_thermal_generation_constructors.jl index 1aa73ee2af..cdb3a38813 100644 --- a/test/test_device_thermal_generation_constructors.jl +++ b/test/test_device_thermal_generation_constructors.jl @@ -906,7 +906,6 @@ end ) solve!(model; output_dir = mktempdir()) - serialize_optimization_model(model, "branch.json") ptdf_vars = get_variable_values(OptimizationProblemResults(model)) power = ptdf_vars[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}( From 25ac4a5f11ad476bcff44c08bed21cf34823a45a Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 15:12:47 -0600 Subject: [PATCH 38/68] fix json3 changes --- src/utils/file_utils.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/utils/file_utils.jl b/src/utils/file_utils.jl index 8ee977d894..ce219a4886 100644 --- a/src/utils/file_utils.jl +++ b/src/utils/file_utils.jl @@ -3,7 +3,7 @@ Return a decoded JSON file. """ function read_json(filename::AbstractString) open(filename, "r") do io - JSON3.parse(io) + JSON3.read(io) end end @@ -28,7 +28,7 @@ end function read_file_hashes(path) data = open(joinpath(path, IS.HASH_FILENAME), "r") do io - JSON3.parse(io) + JSON3.read(io) end return data["files"] From 81727596365ebc3cc15c00cc3cf58f22976f3caf Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 15:12:56 -0600 Subject: [PATCH 39/68] add additional aqua tests --- test/runtests.jl | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index f6eaa26fe9..9c0568af49 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,9 +2,11 @@ include("includes.jl") # Code Quality Tests import Aqua -Aqua.test_unbound_args(PowerSimulations) Aqua.test_undefined_exports(PowerSimulations) Aqua.test_ambiguities(PowerSimulations) +Aqua.test_stale_deps(PowerSimulations) +Aqua.test_persistent_tasks(PowerSimulations) +Aqua.test_unbound_args(PowerSimulations) const LOG_FILE = "power-simulations-test.log" From f0fade0b4d7ebb408960d1e6b96e139fc7140b23 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 15:28:52 -0600 Subject: [PATCH 40/68] update function signature update initial conditions --- ...itial_conditions_update_in_memory_store.jl | 144 +++++++++++++----- .../initial_condition_update_simulation.jl | 114 +++++++++++--- 2 files changed, 198 insertions(+), 60 deletions(-) diff --git a/src/operation/initial_conditions_update_in_memory_store.jl b/src/operation/initial_conditions_update_in_memory_store.jl index 48d56e9afd..36707f3b41 100644 --- a/src/operation/initial_conditions_update_in_memory_store.jl +++ b/src/operation/initial_conditions_update_in_memory_store.jl @@ -2,15 +2,25 @@ ################## ic updates from store for emulation problems simulation ################# function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{InitialTimeDurationOn, Nothing}, - InitialCondition{InitialTimeDurationOn, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialTimeDurationOn, Nothing}, + InitialCondition{InitialTimeDurationOn, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialTimeDurationOn, Nothing}, + InitialCondition{InitialTimeDurationOn, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, TimeDurationOn(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -19,15 +29,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{InitialTimeDurationOff, Nothing}, - InitialCondition{InitialTimeDurationOff, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialTimeDurationOff, Nothing}, + InitialCondition{InitialTimeDurationOff, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialTimeDurationOff, Nothing}, + InitialCondition{InitialTimeDurationOff, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, TimeDurationOff(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -36,15 +56,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{DevicePower, Nothing}, - InitialCondition{DevicePower, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{DevicePower, Nothing}, + InitialCondition{DevicePower, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DevicePower, Nothing}, + InitialCondition{DevicePower, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, ActivePowerVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -53,15 +83,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{DeviceStatus, Nothing}, - InitialCondition{DeviceStatus, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{DeviceStatus, Nothing}, + InitialCondition{DeviceStatus, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DeviceStatus, Nothing}, + InitialCondition{DeviceStatus, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, OnVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) @@ -70,15 +110,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{DeviceAboveMinPower, Nothing}, - InitialCondition{DeviceAboveMinPower, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{DeviceAboveMinPower, Nothing}, + InitialCondition{DeviceAboveMinPower, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DeviceAboveMinPower, Nothing}, + InitialCondition{DeviceAboveMinPower, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, PowerAboveMinimumVariable(), get_component_type(ic)) @@ -104,15 +154,25 @@ end =# function update_initial_conditions!( - ics::Vector{ - Union{ - InitialCondition{InitialEnergyLevel, Nothing}, - InitialCondition{InitialEnergyLevel, JuMP.VariableRef}, - }, - }, + ics::T, store::EmulationModelStore, ::Dates.Millisecond, -) +) where { + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialEnergyLevel, Nothing}, + InitialCondition{InitialEnergyLevel, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialEnergyLevel, Nothing}, + InitialCondition{InitialEnergyLevel, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_value(store, EnergyVariable(), get_component_type(ic)) set_ic_quantity!(ic, get_last_recorded_value(var_val)[get_component_name(ic)]) diff --git a/src/simulation/initial_condition_update_simulation.jl b/src/simulation/initial_condition_update_simulation.jl index 4f4a46559e..48f6c4fd74 100644 --- a/src/simulation/initial_condition_update_simulation.jl +++ b/src/simulation/initial_condition_update_simulation.jl @@ -22,12 +22,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, model_resolution::Dates.Millisecond, ) where { - T <: InitialCondition{InitialTimeDurationOn, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialTimeDurationOn, Nothing}, + InitialCondition{InitialTimeDurationOn, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialTimeDurationOn, Nothing}, + InitialCondition{InitialTimeDurationOn, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_system_state_value(state, TimeDurationOn(), get_component_type(ic)) state_resolution = get_data_resolution( @@ -42,12 +55,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, model_resolution::Dates.Millisecond, ) where { - T <: InitialCondition{InitialTimeDurationOff, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialTimeDurationOff, Nothing}, + InitialCondition{InitialTimeDurationOff, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialTimeDurationOff, Nothing}, + InitialCondition{InitialTimeDurationOff, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_system_state_value(state, TimeDurationOff(), get_component_type(ic)) state_resolution = get_data_resolution( @@ -62,12 +88,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, ::Dates.Millisecond, ) where { - T <: InitialCondition{DevicePower, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{DevicePower, Nothing}, + InitialCondition{DevicePower, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DevicePower, Nothing}, + InitialCondition{DevicePower, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics comp_name = get_component_name(ic) comp_type = get_component_type(ic) @@ -101,12 +140,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, ::Dates.Millisecond, ) where { - T <: InitialCondition{DeviceStatus, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{DeviceStatus, Nothing}, + InitialCondition{DeviceStatus, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DeviceStatus, Nothing}, + InitialCondition{DeviceStatus, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_system_state_value(state, OnVariable(), get_component_type(ic)) set_ic_quantity!(ic, var_val[get_component_name(ic)]) @@ -115,12 +167,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, ::Dates.Millisecond, ) where { - T <: InitialCondition{DeviceAboveMinPower, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{DeviceAboveMinPower, Nothing}, + InitialCondition{DeviceAboveMinPower, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{DeviceAboveMinPower, Nothing}, + InitialCondition{DeviceAboveMinPower, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_system_state_value( state, @@ -133,12 +198,25 @@ function update_initial_conditions!( end function update_initial_conditions!( - ics::Vector{T}, + ics::T, state::SimulationState, ::Dates.Millisecond, ) where { - T <: InitialCondition{InitialEnergyLevel, S}, -} where {S <: Union{Float64, JuMP.VariableRef}} + T <: Union{ + Vector{ + Union{ + InitialCondition{InitialEnergyLevel, Nothing}, + InitialCondition{InitialEnergyLevel, Float64}, + }, + }, + Vector{ + Union{ + InitialCondition{InitialEnergyLevel, Nothing}, + InitialCondition{InitialEnergyLevel, JuMP.VariableRef}, + }, + }, + }, +} for ic in ics var_val = get_system_state_value(state, EnergyVariable(), get_component_type(ic)) set_ic_quantity!(ic, var_val[get_component_name(ic)]) From cfad9decc76e99b1a3943ec7929ad942358e44a4 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 16:10:54 -0600 Subject: [PATCH 41/68] add cross package test --- .github/workflows/cross-package-test.yml | 47 ++++++++++++++++++++++++ 1 file changed, 47 insertions(+) create mode 100644 .github/workflows/cross-package-test.yml diff --git a/.github/workflows/cross-package-test.yml b/.github/workflows/cross-package-test.yml new file mode 100644 index 0000000000..cda8a3ede3 --- /dev/null +++ b/.github/workflows/cross-package-test.yml @@ -0,0 +1,47 @@ +name: CrossPackageTest + +on: + push: + branches: [main] + tags: [v*] + pull_request: + +jobs: + test: + name: Julia v${{ matrix.julia-version }} - ${{ matrix.package_name }} + runs-on: ${{ matrix.os }} + strategy: + matrix: + julia-version: [1] + os: [ubuntu-latest] + package_name: [HydroPowerSimulations, StorageSystemsSimulations] + continue-on-error: true + steps: + - uses: actions/checkout@v2 + - uses: julia-actions/setup-julia@v1 + with: + version: ${{ matrix.julia-version }} + arch: x64 + - uses: julia-actions/julia-buildpkg@latest + - name: Clone ${{matrix.package_name}} + uses: actions/checkout@v2 + with: + repository: NREL-Sienna/${{matrix.package_name}}.jl + path: downstream + - name: Run the tests + shell: julia --project=downstream {0} + run: | + using Pkg + try + # Force it to use this PR's version of the package + Pkg.develop(PackageSpec(path=".")) # resolver may fail with main deps + Pkg.update() + Pkg.test() # resolver may fail with test time deps + catch err + err isa Pkg.Resolve.ResolverError || rethrow() + # If we can't resolve that means this is incompatible by SemVer, and this is fine. + # It means we marked this as a breaking change, so we don't need to worry about + # mistakenly introducing a breaking change as we have intentionally made one. + @info "Not compatible with this release. No problem." exception=err + exit(0) # Exit immediately, as a success + end From ab0de86a9c520a72120e63e01218462cef01bb89 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Wed, 25 Sep 2024 17:19:51 -0600 Subject: [PATCH 42/68] remove debug --- 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 e092cea01b..bc89287be1 100644 --- a/src/operation/operation_model_interface.jl +++ b/src/operation/operation_model_interface.jl @@ -114,7 +114,7 @@ function solve_impl!(model::OperationModel) tss = replace("$(ts)", ":" => "_") model_export_path = joinpath(model_output_dir, "exported_$(model_name)_$(tss).json") serialize_optimization_model(container, model_export_path) - @debug write_lp_file( + write_lp_file( get_jump_model(container), replace(model_export_path, ".json" => ".lp"), ) From 88dfd7f8f2cda4c787507ca24394360f7c4fabb1 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 27 Sep 2024 13:17:46 -0600 Subject: [PATCH 43/68] update tests --- test/test_simulation_execute.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_simulation_execute.jl b/test/test_simulation_execute.jl index 22b8cbb251..d22a2efaba 100644 --- a/test/test_simulation_execute.jl +++ b/test/test_simulation_execute.jl @@ -42,7 +42,7 @@ end test_path = joinpath(folder, "consecutive", "problems", "ED", "optimization_model_exports") @test ispath(test_path) - @test length(readdir(test_path)) == 2 + @test length(readdir(test_path)) == 4 end function test_2_stage_decision_models_with_feedforwards(in_memory) From a3a3f256243b853ab6c35a54b3b6770667c6aab0 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Sat, 28 Sep 2024 16:34:10 -0600 Subject: [PATCH 44/68] fix incorrect no load cost --- .../devices/common/objective_function/common.jl | 3 ++- src/devices_models/devices/thermal_generation.jl | 6 ++++-- 2 files changed, 6 insertions(+), 3 deletions(-) diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index 3d665d35e4..9a8b65ce81 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -54,7 +54,8 @@ function add_proportional_cost!( cost_term = proportional_cost(op_cost_data, U(), d, V()) iszero(cost_term) && continue for t in get_time_steps(container) - _add_proportional_term!(container, U(), d, cost_term * multiplier, t) + exp = _add_proportional_term!(container, U(), d, cost_term * multiplier, t) + add_to_expression!(container, ProductionCostExpression, exp, d, t) end end return diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index bfd514ca80..633de9b9be 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -75,8 +75,10 @@ initial_condition_variable(::InitialTimeDurationOff, d::PSY.ThermalGen, ::Abstra ########################Objective Function################################################## # TODO: Decide what is the cost for OnVariable, if fixed or constant term in variable -#proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractThermalFormulation) = no_load_cost(cost, S, T, U) -proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractThermalFormulation) = PSY.get_fixed(cost) +function proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractThermalFormulation) + return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) +end + proportional_cost(cost::PSY.MarketBidCost, ::OnVariable, ::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_no_load_cost(cost) has_multistart_variables(::PSY.ThermalGen, ::AbstractThermalFormulation)=false From f617917022abb12a5a32accecbfa9b947c0fe67d Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Sat, 28 Sep 2024 17:26:57 -0600 Subject: [PATCH 45/68] add monotonicity check --- .../objective_function/quadratic_curve.jl | 34 ++++++++++++++++++- 1 file changed, 33 insertions(+), 1 deletion(-) diff --git a/src/devices_models/devices/common/objective_function/quadratic_curve.jl b/src/devices_models/devices/common/objective_function/quadratic_curve.jl index 4f5a8aea09..686dc0d236 100644 --- a/src/devices_models/devices/common/objective_function/quadratic_curve.jl +++ b/src/devices_models/devices/common/objective_function/quadratic_curve.jl @@ -45,7 +45,15 @@ function _add_quadraticcurve_variable_cost!( proportional_term_per_unit::Vector{Float64}, quadratic_term_per_unit::Vector{Float64}, ) where {T <: VariableType} + lb, ub = PSY.get_active_power_limits(component) for t in get_time_steps(container) + _check_quadratic_monotonicity( + PSY.get_name(component), + quadratic_term_per_unit[t], + proportional_term_per_unit[t], + lb, + ub, + ) _add_quadraticcurve_variable_term_to_model!( container, T(), @@ -66,6 +74,13 @@ function _add_quadraticcurve_variable_cost!( proportional_term_per_unit::Float64, quadratic_term_per_unit::Float64, ) where {T <: VariableType} + lb, ub = PSY.get_active_power_limits(component) + _check_quadratic_monotonicity(PSY.get_name(component), + quadratic_term_per_unit, + proportional_term_per_unit, + lb, + ub, + ) for t in get_time_steps(container) _add_quadraticcurve_variable_term_to_model!( container, @@ -79,6 +94,23 @@ function _add_quadraticcurve_variable_cost!( return end +function _check_quadratic_monotonicity( + name::String, + quad_term::Float64, + linear_term::Float64, + lb::Float64, + ub::Float64, +) + fp_lb = 2 * quad_term * lb + linear_term + fp_ub = 2 * quad_term * ub + linear_term + + if fp_lb < 0 || fp_ub < 0 + @warn "Cost function for component $name is not monotonically increasing in the range [$lb, $ub]. \ + This can lead to unexpected results" + end + return +end + @doc raw""" Adds to the cost function cost terms for sum of variables with common factor to be used for cost expression for optimization_container model. @@ -149,7 +181,7 @@ function _add_variable_cost_to_objective!( } throw( IS.ConflictingInputsError( - "Quadratic Cost Curves are not allowed for Compact formulations", + "Quadratic Cost Curves are not compatible with Compact formulations", ), ) return From dbd26d3036df1d900f2f22f2dbefc18990b41e83 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 30 Sep 2024 11:56:05 -0600 Subject: [PATCH 46/68] add must run to performance test --- test/performance/performance_test.jl | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/test/performance/performance_test.jl b/test/performance/performance_test.jl index d387796ccf..9c97cce33d 100644 --- a/test/performance/performance_test.jl +++ b/test/performance/performance_test.jl @@ -22,6 +22,11 @@ try sys_rts_rt = build_system(PSISystems, "modified_RTS_GMLC_RT_sys") sys_rts_realization = build_system(PSISystems, "modified_RTS_GMLC_realization_sys") + for sys in [sys_rts_da, sys_rts_rt, sys_rts_realization] + g = get_component(ThermalStandard, sys, "121_NUCLEAR_1") + set_must_run!(g, true) + end + for i in 1:2 template_uc = ProblemTemplate( NetworkModel( @@ -32,8 +37,8 @@ try ), ) - set_device_model!(template_uc, ThermalMultiStart, ThermalCompactUnitCommitment) - set_device_model!(template_uc, ThermalStandard, ThermalCompactUnitCommitment) + set_device_model!(template_uc, ThermalMultiStart, ThermalStandardUnitCommitment) + set_device_model!(template_uc, ThermalStandard, ThermalStandardUnitCommitment) set_device_model!(template_uc, RenewableDispatch, RenewableFullDispatch) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) set_device_model!(template_uc, DeviceModel(Line, StaticBranch)) From fb967c5eaf4f20dd4e64e4e6efe366279c4cdd28 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Mon, 30 Sep 2024 20:10:39 -0600 Subject: [PATCH 47/68] fix pwl implementation for must run --- .../objective_function/piecewise_linear.jl | 67 +++++++++++++------ 1 file changed, 47 insertions(+), 20 deletions(-) diff --git a/src/devices_models/devices/common/objective_function/piecewise_linear.jl b/src/devices_models/devices/common/objective_function/piecewise_linear.jl index bc184fb972..d173cd85d8 100644 --- a/src/devices_models/devices/common/objective_function/piecewise_linear.jl +++ b/src/devices_models/devices/common/objective_function/piecewise_linear.jl @@ -54,6 +54,52 @@ end ################# PWL Constraints ################ ################################################## +function _determine_bin_lhs( + container::OptimizationContainer, + sos_status::SOSStatusVariable, + component::T, + period::Int) where {T <: PSY.Component} + name = PSY.get_name(component) + if sos_status == SOSStatusVariable.NO_VARIABLE + return 1.0 + @debug "Using Piecewise Linear cost function but no variable/parameter ref for ON status is passed. Default status will be set to online (1.0)" _group = + LOG_GROUP_COST_FUNCTIONS + + elseif sos_status == SOSStatusVariable.PARAMETER + param = get_default_on_parameter(component) + return get_parameter(container, param, T).parameter_array[name, period] + @debug "Using Piecewise Linear cost function with parameter OnStatusParameter, $T" _group = + LOG_GROUP_COST_FUNCTIONS + elseif sos_status == SOSStatusVariable.VARIABLE + var = get_default_on_variable(component) + return get_variable(container, var, T)[name, period] + @debug "Using Piecewise Linear cost function with variable OnVariable $T" _group = + LOG_GROUP_COST_FUNCTIONS + else + @assert false + end +end + +function _get_bin_lhs( + container::OptimizationContainer, + sos_status::SOSStatusVariable, + component::T, + period::Int) where {T <: PSY.Component} + return _determine_bin_lhs(container, sos_status, component, period) +end + +function _get_bin_lhs( + container::OptimizationContainer, + sos_status::SOSStatusVariable, + component::PSY.ThermalGen, + period::Int) + if PSY.get_must_run(component) + return 1.0 + else + return _determine_bin_lhs(container, sos_status, component, period) + end +end + """ Implement the constraints for PWL variables. That is: @@ -86,26 +132,7 @@ function _add_pwl_constraint!( variables[name, period] == sum(pwl_vars[name, ix, period] * break_points[ix] for ix in 1:len_cost_data) ) - - if sos_status == SOSStatusVariable.NO_VARIABLE - bin = 1.0 - @debug "Using Piecewise Linear cost function but no variable/parameter ref for ON status is passed. Default status will be set to online (1.0)" _group = - LOG_GROUP_COST_FUNCTIONS - - elseif sos_status == SOSStatusVariable.PARAMETER - param = get_default_on_parameter(component) - bin = get_parameter(container, param, T).parameter_array[name, period] - @debug "Using Piecewise Linear cost function with parameter OnStatusParameter, $T" _group = - LOG_GROUP_COST_FUNCTIONS - elseif sos_status == SOSStatusVariable.VARIABLE - var = get_default_on_variable(component) - bin = get_variable(container, var, T)[name, period] - @debug "Using Piecewise Linear cost function with variable OnVariable $T" _group = - LOG_GROUP_COST_FUNCTIONS - else - @assert false - end - + bin = _get_bin_lhs(container, sos_status, component, period) const_normalization_container = lazy_container_addition!( container, PieceWiseLinearCostConstraint(), From 307497236202724b72b65de06aa8d4204870451e Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 1 Oct 2024 13:53:23 -0700 Subject: [PATCH 48/68] remove SOS formulation --- src/core/formulations.jl | 4 -- .../device_constructors/branch_constructor.jl | 4 +- src/devices_models/devices/AC_branches.jl | 48 ---------------- .../devices/TwoTerminalDC_branches.jl | 56 +------------------ .../devices/common/add_to_expression.jl | 4 +- 5 files changed, 6 insertions(+), 110 deletions(-) diff --git a/src/core/formulations.jl b/src/core/formulations.jl index 49bd495d71..f7ced88017 100644 --- a/src/core/formulations.jl +++ b/src/core/formulations.jl @@ -138,10 +138,6 @@ struct HVDCTwoTerminalDispatch <: AbstractTwoTerminalDCLineFormulation end Branch type to represent piecewise lossy power flow on two terminal DC lines """ struct HVDCTwoTerminalPiecewiseLoss <: AbstractTwoTerminalDCLineFormulation end -""" -Branch type to represent piecewise lossy power flow on two terminal DC lines using SOS2 -""" -struct HVDCTwoTerminalSOSPiecewiseLoss <: AbstractTwoTerminalDCLineFormulation end # Not Implemented # struct VoltageSourceDC <: AbstractTwoTerminalDCLineFormulation end diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 94c6f4bc61..4fb90c8bb3 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -765,7 +765,7 @@ function construct_device!( network_model::NetworkModel{<:AbstractPTDFModel}, ) where { T <: TwoTerminalHVDCTypes, - U <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + U <: HVDCTwoTerminalPiecewiseLoss, } devices = get_available_components(model, sys) @@ -809,7 +809,7 @@ function construct_device!( network_model::NetworkModel{<:AbstractPTDFModel}, ) where { T <: TwoTerminalHVDCTypes, - U <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + U <: HVDCTwoTerminalPiecewiseLoss, } devices = get_available_components(model, sys) diff --git a/src/devices_models/devices/AC_branches.jl b/src/devices_models/devices/AC_branches.jl index e20d1edbfc..8046f17242 100644 --- a/src/devices_models/devices/AC_branches.jl +++ b/src/devices_models/devices/AC_branches.jl @@ -251,54 +251,6 @@ function _add_dense_pwl_loss_variables!( end end -#SOS Model -function _add_sparse_pwl_loss_variables!( - container::OptimizationContainer, - devices, - ::DeviceModel{D, HVDCTwoTerminalSOSPiecewiseLoss}, -) where {D <: TwoTerminalHVDCTypes} - # Check if type and length of PWL loss model are the same for all devices - #_check_pwl_loss_model(devices) - - # Create Variables - time_steps = get_time_steps(container) - settings = get_settings(container) - formulation = HVDCTwoTerminalPiecewiseLoss() - T = HVDCPiecewiseLossVariable - binary = get_variable_binary(T(), D, formulation) - first_loss = PSY.get_loss(first(devices)) - if isa(first_loss, PSY.LinearCurve) - len_segments = 4 # 2*1 + 2 - elseif isa(first_loss, PSY.PiecewiseIncrementalCurve) - len_segments = 2 * length(PSY.get_slopes(first_loss)) + 2 - else - error("Should not be here") - end - - T = HVDCPiecewiseLossVariable - var_container = lazy_container_addition!(container, T(), D) - - for d in devices - name = PSY.get_name(d) - for t in time_steps - pwlvars = Array{JuMP.VariableRef}(undef, len_segments) - for i in 1:len_segments - pwlvars[i] = - var_container[(name, i, t)] = JuMP.@variable( - get_jump_model(container), - base_name = "$(T)_$(name)_{pwl_$(i), $(t)}", - binary = binary - ) - ub = get_variable_upper_bound(T(), d, formulation) - ub !== nothing && JuMP.set_upper_bound(var_container[name, i, t], ub) - - lb = get_variable_lower_bound(T(), d, formulation) - lb !== nothing && JuMP.set_lower_bound(var_container[name, i, t], lb) - end - end - end -end - # Full Binary function _add_sparse_pwl_loss_variables!( container::OptimizationContainer, diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 26acd94857..93864228a3 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -340,58 +340,6 @@ function add_constraints!( return end -# SOS Model -function add_constraints!( - container::OptimizationContainer, - ::Type{T}, - devices::Union{Vector{U}, IS.FlattenIteratorWrapper{U}}, - ::DeviceModel{U, HVDCTwoTerminalSOSPiecewiseLoss}, - ::NetworkModel{<:PM.AbstractPowerModel}, -) where {T <: HVDCFlowCalculationConstraint, U <: PSY.TwoTerminalHVDCLine} - var_pwl = get_variable(container, HVDCPiecewiseLossVariable(), U) - names = PSY.get_name.(devices) - time_steps = get_time_steps(container) - flow_ft = get_variable(container, HVDCActivePowerReceivedFromVariable(), U) - flow_tf = get_variable(container, HVDCActivePowerReceivedToVariable(), U) - - constraint_from_to = - add_constraints_container!(container, T(), U, names, time_steps; meta = "ft") - constraint_to_from = - add_constraints_container!(container, T(), U, names, time_steps; meta = "tf") - for d in devices - name = PSY.get_name(d) - loss = PSY.get_loss(d) - from_to_params, to_from_params = _get_pwl_loss_params(d, loss) - #@show from_to_params - #@show to_from_params - segments = _get_range_segments(d, loss) - for t in time_steps - ## Add Equality Constraints ## - constraint_from_to[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - flow_ft[name, t] == sum( - var_pwl[name, ix, t] * from_to_params[ix] for - ix in segments - ) - ) - constraint_to_from[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - flow_tf[name, t] == sum( - var_pwl[name, ix, t] * to_from_params[ix] for - ix in segments - ) - ) - ## Add SOS Constraints ### - pwl_vars_subset = [var_pwl[name, s, t] for s in segments] - JuMP.@constraint( - get_jump_model(container), - pwl_vars_subset in MOI.SOS2(collect(segments)) - ) - end - end - return -end - #################################### Rate Limits Constraints ################################################## function _get_flow_bounds(d::PSY.TwoTerminalHVDCLine) check_hvdc_line_limits_consistency(d) @@ -594,7 +542,7 @@ function add_constraints!( ) where { T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + V <: HVDCTwoTerminalPiecewiseLoss, } inter_network_branches = U[] for d in devices @@ -633,7 +581,7 @@ function add_constraints!( ) where { T <: Union{FlowRateConstraintFromTo, FlowRateConstraintToFrom}, U <: PSY.TwoTerminalHVDCLine, - V <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + V <: HVDCTwoTerminalPiecewiseLoss, } if T <: FlowRateConstraintFromTo _add_hvdc_flow_constraints!( diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 48d77cc8af..7286427fda 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -389,7 +389,7 @@ function add_to_expression!( T <: ActivePowerBalance, U <: HVDCActivePowerReceivedFromVariable, V <: TwoTerminalHVDCTypes, - W <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + W <: HVDCTwoTerminalPiecewiseLoss, X <: AbstractPTDFModel, } var = get_variable(container, U(), V) @@ -426,7 +426,7 @@ function add_to_expression!( T <: ActivePowerBalance, U <: HVDCActivePowerReceivedToVariable, V <: TwoTerminalHVDCTypes, - W <: Union{HVDCTwoTerminalPiecewiseLoss, HVDCTwoTerminalSOSPiecewiseLoss}, + W <: HVDCTwoTerminalPiecewiseLoss, X <: AbstractPTDFModel, } var = get_variable(container, U(), V) From 7b42fe1c48f94086348caab7dda8ffecc75b27a9 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 1 Oct 2024 15:53:18 -0700 Subject: [PATCH 49/68] update twoterminal dispatch model --- src/core/constraints.jl | 15 --- .../device_constructors/branch_constructor.jl | 1 - .../devices/TwoTerminalDC_branches.jl | 109 ++++++++---------- .../devices/common/add_to_expression.jl | 8 +- 4 files changed, 51 insertions(+), 82 deletions(-) diff --git a/src/core/constraints.jl b/src/core/constraints.jl index a8300d356d..28afddd5cf 100644 --- a/src/core/constraints.jl +++ b/src/core/constraints.jl @@ -332,21 +332,6 @@ The specified constraint is formulated as: ``` """ struct PhaseAngleControlLimit <: ConstraintType end -""" -Struct to create the constraints that set the losses through a lossy HVDC two-terminal line. - -For more information check [Branch Formulations](@ref PowerSystems.Branch-Formulations). - -The specified constraints are formulated as: - -```math -\\begin{align*} -& f_t^\\text{to-from} - f_t^\\text{from-to} \\le \\ell_t,\\quad \\forall t \\in \\{1,\\dots, T\\} \\\\ -& f_t^\\text{from-to} - f_t^\\text{to-from} \\le \\ell_t,\\quad \\forall t \\in \\{1,\\dots, T\\} -\\end{align*} -``` -""" -struct HVDCLossesAbsoluteValue <: ConstraintType end struct HVDCDirection <: ConstraintType end struct InterfaceFlowLimit <: ConstraintType end struct HVDCFlowCalculationConstraint <: ConstraintType end diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index 4fb90c8bb3..b1926fefbb 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -695,7 +695,6 @@ function construct_device!( add_constraints!(container, FlowRateConstraintFromTo, devices, model, network_model) add_constraints!(container, FlowRateConstraintToFrom, devices, model, network_model) add_constraints!(container, HVDCPowerBalance, devices, model, network_model) - add_constraints!(container, HVDCLossesAbsoluteValue, devices, model, network_model) return end diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 93864228a3..b2a0372c14 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -661,6 +661,7 @@ function add_constraints!( tf_var = get_variable(container, FlowActivePowerToFromVariable(), T) ft_var = get_variable(container, FlowActivePowerFromToVariable(), T) direction_var = get_variable(container, HVDCFlowDirectionVariable(), T) + losses = get_variable(container, HVDCLosses(), T) constraint_ft_ub = add_constraints_container!( container, @@ -694,6 +695,30 @@ function add_constraints!( time_steps; meta = "ft_lb", ) + constraint_loss = add_constraints_container!( + container, + HVDCPowerBalance(), + T, + names, + time_steps; + meta = "loss", + ) + constraint_loss_aux1 = add_constraints_container!( + container, + HVDCPowerBalance(), + T, + names, + time_steps; + meta = "loss_aux1", + ) + constraint_loss_aux2 = add_constraints_container!( + container, + HVDCPowerBalance(), + T, + names, + time_steps; + meta = "loss_aux2", + ) for d in devices loss = PSY.get_loss(d) if !isa(loss, PSY.LinearCurve) @@ -704,78 +729,38 @@ function add_constraints!( l1 = PSY.get_proportional_term(loss) l0 = PSY.get_constant_term(loss) name = PSY.get_name(d) + R_min_from, R_max_from = PSY.get_active_power_limits_from(d) + R_min_to, R_max_to = PSY.get_active_power_limits_to(d) for t in get_time_steps(container) - if l1 == 0.0 && l0 == 0.0 - constraint_tf_ub[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - tf_var[name, t] - ft_var[name, t] == 0.0 - ) - constraint_ft_ub[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - ft_var[name, t] - tf_var[name, t] == 0.0 - ) - else - constraint_tf_ub[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - tf_var[name, t] - ft_var[name, t] <= l1 * tf_var[name, t] - l0 - ) - constraint_ft_ub[PSY.get_name(d), t] = JuMP.@constraint( - get_jump_model(container), - ft_var[name, t] - tf_var[name, t] >= l1 * ft_var[name, t] + l0 - ) - end + constraint_tf_ub[PSY.get_name(d), t] = JuMP.@constraint( + get_jump_model(container), + tf_var[name, t] <= R_max_to * direction_var[name, t] + ) constraint_tf_lb[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), - ft_var[name, t] - tf_var[name, t] >= - -M_VALUE * (1 - direction_var[name, t]) + tf_var[name, t] >= R_min_to * (1 - direction_var[name, t]) + ) + constraint_ft_ub[PSY.get_name(d), t] = JuMP.@constraint( + get_jump_model(container), + ft_var[name, t] <= R_max_from * (1 - direction_var[name, t]) ) constraint_ft_lb[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), - tf_var[name, t] - ft_var[name, t] >= -M_VALUE * (direction_var[name, t]) + ft_var[name, t] >= R_min_from * direction_var[name, t] ) - end - end - return -end - -function add_constraints!( - container::OptimizationContainer, - ::Type{HVDCLossesAbsoluteValue}, - devices::IS.FlattenIteratorWrapper{T}, - ::DeviceModel{T, <:AbstractTwoTerminalDCLineFormulation}, - ::NetworkModel{<:PM.AbstractDCPModel}, -) where {T <: PSY.TwoTerminalHVDCLine} - time_steps = get_time_steps(container) - names = [PSY.get_name(d) for d in devices] - losses = get_variable(container, HVDCLosses(), T) - tf_var = get_variable(container, FlowActivePowerToFromVariable(), T) - ft_var = get_variable(container, FlowActivePowerFromToVariable(), T) - constraint_tf = add_constraints_container!( - container, - HVDCLossesAbsoluteValue(), - T, - names, - time_steps; - meta = "tf", - ) - constraint_ft = add_constraints_container!( - container, - HVDCLossesAbsoluteValue(), - T, - names, - time_steps; - meta = "ft", - ) - for d in devices - name = PSY.get_name(d) - for t in get_time_steps(container) - constraint_tf[PSY.get_name(d), t] = JuMP.@constraint( + constraint_loss[PSY.get_name(d), t] = JuMP.@constraint( + get_jump_model(container), + tf_var[name, t] + ft_var[name, t] == losses[name, t] + ) + constraint_loss_aux1[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), - tf_var[name, t] - ft_var[name, t] <= losses[name, t] + losses[name, t] >= + l1 * ft_var[name, t] + l0 - M_VALUE * direction_var[name, t] ) - constraint_ft[PSY.get_name(d), t] = JuMP.@constraint( + constraint_loss_aux2[PSY.get_name(d), t] = JuMP.@constraint( get_jump_model(container), - -tf_var[name, t] + ft_var[name, t] <= losses[name, t] + losses[name, t] >= + l1 * tf_var[name, t] + l0 - M_VALUE * (1 - direction_var[name, t]) ) end end diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 7286427fda..608d28afee 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -329,9 +329,9 @@ function add_to_expression!( ref_bus_to = get_reference_bus(network_model, PSY.get_arc(d).to) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] - _add_to_jump_expression!(nodal_expr[bus_no_to, t], flow_variable, 1.0) + _add_to_jump_expression!(nodal_expr[bus_no_to, t], flow_variable, -1.0) if ref_bus_from != ref_bus_to - _add_to_jump_expression!(sys_expr[ref_bus_to, t], flow_variable, 1.0) + _add_to_jump_expression!(sys_expr[ref_bus_to, t], flow_variable, -1.0) end end end @@ -366,9 +366,9 @@ function add_to_expression!( ref_bus_from = get_reference_bus(network_model, PSY.get_arc(d).from) for t in get_time_steps(container) flow_variable = var[PSY.get_name(d), t] - _add_to_jump_expression!(nodal_expr[bus_no_from, t], flow_variable, 1.0) + _add_to_jump_expression!(nodal_expr[bus_no_from, t], flow_variable, -1.0) if ref_bus_from != ref_bus_to - _add_to_jump_expression!(sys_expr[ref_bus_from, t], flow_variable, 1.0) + _add_to_jump_expression!(sys_expr[ref_bus_from, t], flow_variable, -1.0) end end end From 2fa41a2be3c6c4497674cee2872fdde964b5fa97 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 1 Oct 2024 16:47:54 -0700 Subject: [PATCH 50/68] remove HVDC losses in different asynch regions --- .../devices/common/add_to_expression.jl | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index 608d28afee..51c33b245b 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -260,13 +260,17 @@ function add_to_expression!( for d in devices name = PSY.get_name(d) device_bus_from = PSY.get_arc(d).from + device_bus_to = PSY.get_arc(d).to ref_bus_from = get_reference_bus(network_model, device_bus_from) - for t in get_time_steps(container) - _add_to_jump_expression!( - expression[ref_bus_from, t], - variable[name, t], - get_variable_multiplier(U(), d, W()), - ) + ref_bus_to = get_reference_bus(network_model, device_bus_to) + if ref_bus_from == ref_bus_to + for t in get_time_steps(container) + _add_to_jump_expression!( + expression[ref_bus_from, t], + variable[name, t], + get_variable_multiplier(U(), d, W()), + ) + end end end return From de63f6fc24ebbac0874dbd4e4db5936143d31a9c Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 1 Oct 2024 16:48:59 -0700 Subject: [PATCH 51/68] bump PSY version --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 01523f9175..5a0824e987 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,7 @@ Logging = "1" MathOptInterface = "1" PowerModels = "^0.21" PowerNetworkMatrices = "^0.11" -PowerSystems = "4" +PowerSystems = "^4.4" PrettyTables = "2" ProgressMeter = "^1.5" Serialization = "1" From 4d3770ed8e1494da37507ba3c1531527152fd967 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 2 Oct 2024 16:53:30 -0700 Subject: [PATCH 52/68] update PM translator --- src/network_models/pm_translator.jl | 12 ++++++------ test/test_device_branch_constructors.jl | 4 ++-- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/src/network_models/pm_translator.jl b/src/network_models/pm_translator.jl index 58381f732e..171297715e 100644 --- a/src/network_models/pm_translator.jl +++ b/src/network_models/pm_translator.jl @@ -271,14 +271,14 @@ function get_branch_to_pm( ::Type{<:PM.AbstractDCPModel}, ) PM_branch = Dict{String, Any}( - "loss1" => PSY.get_loss(branch).l1, + "loss1" => PSY.get_proportional_term(PSY.get_loss(branch)), "mp_pmax" => PSY.get_reactive_power_limits_from(branch).max, "model" => 2, "shutdown" => 0.0, "pmaxt" => PSY.get_active_power_limits_to(branch).max, "pmaxf" => PSY.get_active_power_limits_from(branch).max, "startup" => 0.0, - "loss0" => PSY.get_loss(branch).l0, + "loss0" => PSY.get_constant_term(PSY.get_loss(branch)), "pt" => 0.0, "vt" => PSY.get_magnitude(PSY.get_arc(branch).to), "qmaxf" => PSY.get_reactive_power_limits_from(branch).max, @@ -310,14 +310,14 @@ function get_branch_to_pm( ) check_hvdc_line_limits_unidirectional(branch) PM_branch = Dict{String, Any}( - "loss1" => PSY.get_loss(branch).l1, + "loss1" => PSY.get_proportional_term(PSY.get_loss(branch)), "mp_pmax" => PSY.get_reactive_power_limits_from(branch).max, "model" => 2, "shutdown" => 0.0, "pmaxt" => PSY.get_active_power_limits_to(branch).max, "pmaxf" => PSY.get_active_power_limits_from(branch).max, "startup" => 0.0, - "loss0" => PSY.get_loss(branch).l0, + "loss0" => PSY.get_constant_term(PSY.get_loss(branch)), "pt" => 0.0, "vt" => PSY.get_magnitude(PSY.get_arc(branch).to), "qmaxf" => PSY.get_reactive_power_limits_from(branch).max, @@ -348,14 +348,14 @@ function get_branch_to_pm( ::Type{<:PM.AbstractPowerModel}, ) PM_branch = Dict{String, Any}( - "loss1" => PSY.get_loss(branch).l1, + "loss1" => PSY.get_proportional_term(PSY.get_loss(branch)), "mp_pmax" => PSY.get_reactive_power_limits_from(branch).max, "model" => 2, "shutdown" => 0.0, "pmaxt" => PSY.get_active_power_limits_to(branch).max, "pmaxf" => PSY.get_active_power_limits_from(branch).max, "startup" => 0.0, - "loss0" => PSY.get_loss(branch).l0, + "loss0" => PSY.get_constant_term(PSY.get_loss(branch)), "pt" => 0.0, "vt" => PSY.get_magnitude(PSY.get_arc(branch).to), "qmaxf" => PSY.get_reactive_power_limits_from(branch).max, diff --git a/test/test_device_branch_constructors.jl b/test/test_device_branch_constructors.jl index 118fda9c60..2dc350e14f 100644 --- a/test/test_device_branch_constructors.jl +++ b/test/test_device_branch_constructors.jl @@ -323,7 +323,7 @@ end add_component!(sys_5, hvdc) for net_model in [DCPPowerModel, PTDFPowerModel] @testset "$net_model" begin - PSY.set_loss!(hvdc, (l0 = 0.0, l1 = 0.0)) + PSY.set_loss!(hvdc, PSY.LinearCurve(0.0)) template_uc = ProblemTemplate( NetworkModel(net_model; use_slacks = true), ) @@ -431,7 +431,7 @@ end @test isapprox(no_loss_total_gen, ref_total_gen; atol = 0.1) - PSY.set_loss!(hvdc, (l0 = 0.1, l1 = 0.005)) + PSY.set_loss!(hvdc, PSY.LinearCurve(0.005, 0.1)) model_wl = DecisionModel( template_uc, From 54d1d91753a2d81bc4f151aab75a64e60b66a94e Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 2 Oct 2024 16:53:44 -0700 Subject: [PATCH 53/68] remove Direction Constraint --- src/core/constraints.jl | 1 - .../devices/TwoTerminalDC_branches.jl | 48 ------------------- 2 files changed, 49 deletions(-) diff --git a/src/core/constraints.jl b/src/core/constraints.jl index 28afddd5cf..817130ba45 100644 --- a/src/core/constraints.jl +++ b/src/core/constraints.jl @@ -332,7 +332,6 @@ The specified constraint is formulated as: ``` """ struct PhaseAngleControlLimit <: ConstraintType end -struct HVDCDirection <: ConstraintType end struct InterfaceFlowLimit <: ConstraintType end struct HVDCFlowCalculationConstraint <: ConstraintType end diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index b2a0372c14..53aee2bae9 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -601,54 +601,6 @@ function add_constraints!( return end -function add_constraints!( - container::OptimizationContainer, - ::Type{T}, - devices::IS.FlattenIteratorWrapper{U}, - ::DeviceModel{U, HVDCTwoTerminalDispatch}, - ::NetworkModel{<:PM.AbstractDCPModel}, -) where {T <: HVDCDirection, U <: PSY.TwoTerminalHVDCLine} - time_steps = get_time_steps(container) - names = [PSY.get_name(d) for d in devices] - - tf_var = get_variable(container, FlowActivePowerToFromVariable(), U) - ft_var = get_variable(container, FlowActivePowerFromToVariable(), U) - direction_var = get_variable(container, HVDCFlowDirectionVariable(), U) - - constraint_ft_ub = - add_constraints_container!(container, T(), U, names, time_steps; meta = "ft_ub") - constraint_tf_ub = - add_constraints_container!(container, T(), U, names, time_steps; meta = "tf_ub") - constraint_ft_lb = - add_constraints_container!(container, T(), U, names, time_steps; meta = "ft_lb") - constraint_tf_lb = - add_constraints_container!(container, T(), U, names, time_steps; meta = "tf_lb") - for d in devices - min_rate_to, max_rate_to = PSY.get_active_power_limits_to(d) - min_rate_from, max_rate_from = PSY.get_active_power_limits_to(d) - name = PSY.get_name(d) - for t in time_steps - constraint_tf_ub[name, t] = JuMP.@constraint( - get_jump_model(container), - tf_var[name, t] <= max_rate_to * (1 - direction_var[name, t]) - ) - constraint_ft_ub[name, t] = JuMP.@constraint( - get_jump_model(container), - ft_var[name, t] <= max_rate_from * (1 - direction_var[name, t]) - ) - constraint_tf_lb[name, t] = JuMP.@constraint( - get_jump_model(container), - direction_var[name, t] * min_rate_to <= tf_var[name, t] - ) - constraint_ft_lb[name, t] = JuMP.@constraint( - get_jump_model(container), - direction_var[name, t] * min_rate_from <= tf_var[name, t] - ) - end - end - return -end - function add_constraints!( container::OptimizationContainer, ::Type{HVDCPowerBalance}, From 1ab24845309c48dbdbe89fa1d413352ed034ce14 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Wed, 2 Oct 2024 16:54:04 -0700 Subject: [PATCH 54/68] add HVDC losses to DCP model --- src/devices_models/device_constructors/branch_constructor.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/devices_models/device_constructors/branch_constructor.jl b/src/devices_models/device_constructors/branch_constructor.jl index a4a517ce7d..98d90e510f 100644 --- a/src/devices_models/device_constructors/branch_constructor.jl +++ b/src/devices_models/device_constructors/branch_constructor.jl @@ -718,6 +718,7 @@ function construct_device!( HVDCTwoTerminalDispatch(), ) add_variables!(container, HVDCFlowDirectionVariable, devices, HVDCTwoTerminalDispatch()) + add_variables!(container, HVDCLosses, devices, HVDCTwoTerminalDispatch()) add_to_expression!( container, ActivePowerBalance, @@ -893,7 +894,6 @@ function construct_device!( add_constraints!(container, FlowRateConstraintFromTo, devices, model, network_model) add_constraints!(container, FlowRateConstraintToFrom, devices, model, network_model) add_constraints!(container, HVDCPowerBalance, devices, model, network_model) - add_constraints!(container, HVDCDirection, devices, model, network_model) add_constraint_dual!(container, sys, model) add_feedforward_constraints!(container, model, devices) return From 48b82c63212e0b4da6a7f0e1fefe35ce0cf4aaf3 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:19:29 -0600 Subject: [PATCH 55/68] add missing methods for ics --- src/core/initial_conditions.jl | 10 ++++++++-- src/initial_conditions/calculate_initial_condition.jl | 7 +++++++ 2 files changed, 15 insertions(+), 2 deletions(-) diff --git a/src/core/initial_conditions.jl b/src/core/initial_conditions.jl index 5474275422..8cb845430d 100644 --- a/src/core/initial_conditions.jl +++ b/src/core/initial_conditions.jl @@ -39,16 +39,22 @@ function get_condition( return jump_value(p.value) end +function get_condition( + ::InitialCondition{T, Nothing}, +) where {T <: InitialConditionType} + return nothing +end + get_component(ic::InitialCondition) = ic.component get_value(ic::InitialCondition) = ic.value get_component_name(ic::InitialCondition) = PSY.get_name(ic.component) get_component_type(ic::InitialCondition) = typeof(ic.component) get_ic_type( ::Type{InitialCondition{T, U}}, -) where {T <: InitialConditionType, U <: Union{JuMP.VariableRef, Float64}} = T +) where {T <: InitialConditionType, U <: Union{JuMP.VariableRef, Float64, Nothing}} = T get_ic_type( ::InitialCondition{T, U}, -) where {T <: InitialConditionType, U <: Union{JuMP.VariableRef, Float64}} = T +) where {T <: InitialConditionType, U <: Union{JuMP.VariableRef, Float64, Nothing}} = T """ Stores data to populate initial conditions before the build call diff --git a/src/initial_conditions/calculate_initial_condition.jl b/src/initial_conditions/calculate_initial_condition.jl index 62af678391..be5add50ff 100644 --- a/src/initial_conditions/calculate_initial_condition.jl +++ b/src/initial_conditions/calculate_initial_condition.jl @@ -23,3 +23,10 @@ function set_ic_quantity!( ic.value = var_value return end + +function set_ic_quantity!( + ::InitialCondition{T, Nothing}, + ::Float64, +) where {T <: InitialConditionType} + return +end From f15d14a871bbe6acef6c82127b9f57ea843921f8 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:20:15 -0600 Subject: [PATCH 56/68] handle thermal modeling --- .../devices/thermal_generation.jl | 31 +++++++++++++------ 1 file changed, 21 insertions(+), 10 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index 633de9b9be..ae80274519 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -76,7 +76,11 @@ initial_condition_variable(::InitialTimeDurationOff, d::PSY.ThermalGen, ::Abstra ########################Objective Function################################################## # TODO: Decide what is the cost for OnVariable, if fixed or constant term in variable function proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractThermalFormulation) - return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + PSY.get_fixed(cost) +end + +function proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractCompactUnitCommitment) + return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + PSY.get_fixed(cost) end proportional_cost(cost::PSY.MarketBidCost, ::OnVariable, ::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_no_load_cost(cost) @@ -130,8 +134,10 @@ function _no_load_cost(cost_function::Union{PSY.CostCurve{PSY.LinearCurve}, PSY. end function _no_load_cost(cost_function::PSY.FuelCurve{PSY.PiecewisePointCurve}, d::PSY.ThermalGen) - # value_curve = PSY.get_value_curve(cost_function) - # cost = PSY.get_function_data(value_curve) + return 0.0 +end + +function _no_load_cost(cost_function::PSY.FuelCurve{PSY.PiecewiseIncrementalCurve}, d::PSY.ThermalGen) return 0.0 end @@ -826,8 +832,10 @@ function calculate_aux_variable_value!( if isnothing(get_value(ini_cond[ix])) sum_on_var = time_steps[end] else + on_var_name = get_component_name(ini_cond[ix]) ini_cond_value = get_condition(ini_cond[ix]) - on_var = jump_value.(on_variable_results.data[ix, :]) + # On Var doesn't exist for a unit that has must_run = true + on_var = jump_value.(on_variable_results[on_var_name, :]) aux_variable_container.data[ix, :] .= ini_cond_value sum_on_var = sum(on_var) end @@ -864,15 +872,18 @@ function calculate_aux_variable_value!( time_steps = get_time_steps(container) for ix in eachindex(JuMP.axes(aux_variable_container)[1]) - # if its nothing it means the thermal unit was on must run + # if its nothing it means the thermal unit was on must_run = true # so there is nothing to do but continue if isnothing(get_value(ini_cond[ix])) - continue + sum_on_var = 0.0 + else + on_var_name = get_component_name(ini_cond[ix]) + # On Var doesn't exist for a unit that has must run + on_var = jump_value.(on_variable_results[on_var_name, :]) + ini_cond_value = get_condition(ini_cond[ix]) + aux_variable_container.data[ix, :] .= ini_cond_value + sum_on_var = sum(on_var) end - on_var = jump_value.(on_variable_results.data[ix, :]) - ini_cond_value = get_condition(ini_cond[ix]) - aux_variable_container.data[ix, :] .= ini_cond_value - sum_on_var = sum(on_var) if sum_on_var == time_steps[end] # Unit was always on aux_variable_container.data[ix, :] .= 0.0 elseif sum_on_var == 0.0 # Unit was always off From 33036f7a309a3510e0d714e9ac3a4eb74cd8d4a0 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:20:29 -0600 Subject: [PATCH 57/68] update ramp constraints --- .../common/rateofchange_constraints.jl | 46 ++++++++++++++++--- 1 file changed, 39 insertions(+), 7 deletions(-) diff --git a/src/devices_models/devices/common/rateofchange_constraints.jl b/src/devices_models/devices/common/rateofchange_constraints.jl index b5558b56b1..bca3b5ed5f 100644 --- a/src/devices_models/devices/common/rateofchange_constraints.jl +++ b/src/devices_models/devices/common/rateofchange_constraints.jl @@ -195,14 +195,13 @@ function add_semicontinuous_ramp_constraints!( T::Type{<:ConstraintType}, U::Type{S}, devices::IS.FlattenIteratorWrapper{V}, - model::DeviceModel{V, W}, - X::Type{<:PM.AbstractPowerModel}, + ::DeviceModel{V, W}, + ::Type{<:PM.AbstractPowerModel}, ) where { S <: Union{PowerAboveMinimumVariable, ActivePowerVariable}, V <: PSY.Component, W <: AbstractDeviceFormulation, } - parameters = built_for_recurrent_solves(container) time_steps = get_time_steps(container) variable = get_variable(container, U(), V) varstart = get_variable(container, StartVariable(), V) @@ -222,6 +221,7 @@ function add_semicontinuous_ramp_constraints!( add_constraints_container!(container, T(), V, set_name, time_steps; meta = "dn") for ic in initial_conditions_power + component = get_component(ic) name = get_component_name(ic) # This is to filter out devices that dont need a ramping constraint name ∉ set_name && continue @@ -231,26 +231,58 @@ function add_semicontinuous_ramp_constraints!( ic_power = get_value(ic) @debug "add rate_of_change_constraint" name ic_power + if hasmethod(PSY.get_must_run, Tuple{V}) + must_run = PSY.get_must_run(component) + else + must_run = false + end + + if must_run + rhs_up = ramp_limits.up * minutes_per_period + rhd_dn = ramp_limits.down * minutes_per_period + else + rhs_up = + ramp_limits.up * minutes_per_period + power_limits.min * varstart[name, 1] + rhd_dn = + ramp_limits.down * minutes_per_period + power_limits.min * varstop[name, 1] + end + con_up[name, 1] = JuMP.@constraint( get_jump_model(container), expr_up[name, 1] - ic_power <= - ramp_limits.up * minutes_per_period + power_limits.min * varstart[name, 1] + if must_run + ramp_limits.up * minutes_per_period + else + ramp_limits.up * minutes_per_period + power_limits.min * varstart[name, 1] + end ) con_down[name, 1] = JuMP.@constraint( get_jump_model(container), ic_power - expr_dn[name, 1] <= - ramp_limits.down * minutes_per_period + power_limits.min * varstop[name, 1] + if must_run + ramp_limits.down * minutes_per_period + else + ramp_limits.down * minutes_per_period + power_limits.min * varstop[name, 1] + end ) for t in time_steps[2:end] con_up[name, t] = JuMP.@constraint( get_jump_model(container), expr_up[name, t] - variable[name, t - 1] <= - ramp_limits.up * minutes_per_period + power_limits.min * varstart[name, t] + if must_run + ramp_limits.up * minutes_per_period + else + ramp_limits.up * minutes_per_period + power_limits.min * varstart[name, t] + end ) con_down[name, t] = JuMP.@constraint( get_jump_model(container), variable[name, t - 1] - expr_dn[name, t] <= - ramp_limits.down * minutes_per_period + power_limits.min * varstop[name, t] + if must_run + ramp_limits.down * minutes_per_period + else + ramp_limits.down * minutes_per_period + power_limits.min * varstop[name, t] + end ) end end From d7def69567aa3c7e00cd6ca1f1d2ac83c2f44acb Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:20:53 -0600 Subject: [PATCH 58/68] add must run conditionals --- .../devices/common/add_to_expression.jl | 5 +++- .../common/objective_function/common.jl | 16 +++++++--- src/feedforward/feedforward_constraints.jl | 30 ++++++++++++------- 3 files changed, 36 insertions(+), 15 deletions(-) diff --git a/src/devices_models/devices/common/add_to_expression.jl b/src/devices_models/devices/common/add_to_expression.jl index ca464980ae..e50d55294c 100644 --- a/src/devices_models/devices/common/add_to_expression.jl +++ b/src/devices_models/devices/common/add_to_expression.jl @@ -1269,6 +1269,9 @@ function add_to_expression!( end expression = get_expression(container, T(), V) for d in devices + if PSY.get_must_run(d) + continue + end mult = get_expression_multiplier(U(), T(), d, W()) for t in get_time_steps(container) name = PSY.get_name(d) @@ -1419,7 +1422,7 @@ end function add_to_expression!( container::OptimizationContainer, ::Type{S}, - cost_expression::JuMP.AbstractJuMPScalar, + cost_expression::Union{JuMP.AbstractJuMPScalar, Float64}, component::T, time_period::Int, ) where {S <: CostExpressions, T <: PSY.Component} diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index 9a8b65ce81..c2da8f5a3a 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -41,7 +41,6 @@ end ################################## ####### Proportional Cost ######## ################################## - function add_proportional_cost!( container::OptimizationContainer, ::U, @@ -70,15 +69,24 @@ function add_proportional_cost!( ::U, devices::IS.FlattenIteratorWrapper{T}, ::V, -) where {T <: PSY.ThermalGen, U <: OnVariable, V <: AbstractCompactUnitCommitment} +) where {T <: PSY.ThermalGen, U <: OnVariable, V <: AbstractThermalUnitCommitment} multiplier = objective_function_multiplier(U(), V()) for d in devices op_cost_data = PSY.get_operation_cost(d) cost_term = proportional_cost(op_cost_data, U(), d, V()) iszero(cost_term) && continue for t in get_time_steps(container) - exp = _add_proportional_term!(container, U(), d, cost_term * multiplier, t) - add_to_expression!(container, ProductionCostExpression, exp, d, t) + if !PSY.get_must_run(d) + exp = _add_proportional_term!(container, U(), d, cost_term * multiplier, t) + add_to_expression!(container, ProductionCostExpression, exp, d, t) + end + add_to_expression!( + container, + ProductionCostExpression, + cost_term * multiplier, + d, + t, + ) end end return diff --git a/src/feedforward/feedforward_constraints.jl b/src/feedforward/feedforward_constraints.jl index a85333b85f..656c2deb89 100644 --- a/src/feedforward/feedforward_constraints.jl +++ b/src/feedforward/feedforward_constraints.jl @@ -178,12 +178,17 @@ function _lower_bound_range_with_parameter!( devices::IS.FlattenIteratorWrapper{V}, ) where {V <: PSY.Component} time_steps = axes(constraint_container)[2] - for device in devices, t in time_steps + for device in devices + if hasmethod(PSY.get_must_run, Tuple{V}) + PSY.get_must_run(device) && continue + end name = PSY.get_name(device) - constraint_container[name, t] = JuMP.@constraint( - jump_model, - lhs_array[name, t] >= param_multiplier[name, t] * param_array[name, t] - ) + for t in time_steps + constraint_container[name, t] = JuMP.@constraint( + jump_model, + lhs_array[name, t] >= param_multiplier[name, t] * param_array[name, t] + ) + end end return end @@ -197,12 +202,17 @@ function _upper_bound_range_with_parameter!( devices::IS.FlattenIteratorWrapper{V}, ) where {V <: PSY.Component} time_steps = axes(constraint_container)[2] - for device in devices, t in time_steps + for device in devices + if hasmethod(PSY.get_must_run, Tuple{V}) + PSY.get_must_run(device) && continue + end name = PSY.get_name(device) - constraint_container[name, t] = JuMP.@constraint( - jump_model, - lhs_array[name, t] <= param_multiplier[name, t] * param_array[name, t] - ) + for t in time_steps + constraint_container[name, t] = JuMP.@constraint( + jump_model, + lhs_array[name, t] <= param_multiplier[name, t] * param_array[name, t] + ) + end end return end From 4d537ac6f53cc35818b181f6440867e0fd04e089 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:20:59 -0600 Subject: [PATCH 59/68] fix recorder --- src/utils/recorder_events.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/utils/recorder_events.jl b/src/utils/recorder_events.jl index 755818c522..57d3704141 100644 --- a/src/utils/recorder_events.jl +++ b/src/utils/recorder_events.jl @@ -60,8 +60,8 @@ end function InitialConditionUpdateEvent( simulation_time, ic::InitialCondition, - previous_value::Float64, - model_name, + previous_value::Union{Nothing, Float64}, + model_name::Symbol, ) return InitialConditionUpdateEvent( IS.RecorderEventCommon("InitialConditionUpdateEvent"), @@ -69,8 +69,8 @@ function InitialConditionUpdateEvent( string(get_ic_type(ic)), string(get_component_type(ic)), get_component_name(ic), - get_condition(ic), - previous_value, + isnothing(get_condition(ic)) ? 1e8 : get_condition(ic), + isnothing(previous_value) ? 1e8 : previous_value, string(model_name), ) end From e4cf9ee48e4fc7bb18da43272dafa082bd240de5 Mon Sep 17 00:00:00 2001 From: Jose Daniel Lara Date: Fri, 4 Oct 2024 16:21:10 -0600 Subject: [PATCH 60/68] handle parameters and initial conditions --- .../add_initial_condition.jl | 2 +- src/parameters/add_parameters.jl | 46 +++++++++++++++++++ .../initial_condition_update_simulation.jl | 10 +++- 3 files changed, 55 insertions(+), 3 deletions(-) diff --git a/src/initial_conditions/add_initial_condition.jl b/src/initial_conditions/add_initial_condition.jl index 6c171a7763..a5d69f0483 100644 --- a/src/initial_conditions/add_initial_condition.jl +++ b/src/initial_conditions/add_initial_condition.jl @@ -237,7 +237,7 @@ function add_initial_condition!( ) where { T <: PSY.ThermalGen, U <: AbstractThermalFormulation, - D <: InitialConditionType, + D <: Union{InitialTimeDurationOff, InitialTimeDurationOn, DeviceStatus}, } if get_rebuild_model(get_settings(container)) && has_container_key(container, D, T) return diff --git a/src/parameters/add_parameters.jl b/src/parameters/add_parameters.jl index 6261ded49c..76f1ef2095 100644 --- a/src/parameters/add_parameters.jl +++ b/src/parameters/add_parameters.jl @@ -338,6 +338,52 @@ function _add_parameters!( return end +function _add_parameters!( + container::OptimizationContainer, + ::T, + key::VariableKey{U, D}, + model::DeviceModel{D, W}, + devices::V, +) where { + T <: OnStatusParameter, + U <: OnVariable, + V <: Union{Vector{D}, IS.FlattenIteratorWrapper{D}}, + W <: AbstractThermalFormulation, +} where {D <: PSY.ThermalGen} + @debug "adding" T D U _group = LOG_GROUP_OPTIMIZATION_CONTAINER + names = [PSY.get_name(device) for device in devices if !PSY.get_must_run(device)] + time_steps = get_time_steps(container) + parameter_container = add_param_container!(container, T(), D, key, names, time_steps) + jump_model = get_jump_model(container) + for d in devices + if PSY.get_must_run(d) + continue + end + name = PSY.get_name(d) + if get_variable_warm_start_value(U(), d, W()) === nothing + inital_parameter_value = 0.0 + else + inital_parameter_value = get_variable_warm_start_value(U(), d, W()) + end + for t in time_steps + set_multiplier!( + parameter_container, + get_parameter_multiplier(T(), d, W()), + name, + t, + ) + set_parameter!( + parameter_container, + jump_model, + inital_parameter_value, + name, + t, + ) + end + end + return +end + function _add_parameters!( container::OptimizationContainer, ::T, diff --git a/src/simulation/initial_condition_update_simulation.jl b/src/simulation/initial_condition_update_simulation.jl index 48f6c4fd74..d360753019 100644 --- a/src/simulation/initial_condition_update_simulation.jl +++ b/src/simulation/initial_condition_update_simulation.jl @@ -75,6 +75,7 @@ function update_initial_conditions!( }, } for ic in ics + isnothing(get_value(ic)) && continue var_val = get_system_state_value(state, TimeDurationOff(), get_component_type(ic)) state_resolution = get_data_resolution( get_system_state_data(state, TimeDurationOff(), get_component_type(ic)), @@ -110,10 +111,14 @@ function update_initial_conditions!( for ic in ics comp_name = get_component_name(ic) comp_type = get_component_type(ic) - status_val = get_system_state_value(state, OnVariable(), comp_type)[comp_name] + comp = get_component(ic) + if hasmethod(PSY.get_must_run, Tuple{comp_type}) && PSY.get_must_run(comp) + status_val = 1.0 + else + status_val = get_system_state_value(state, OnVariable(), comp_type)[comp_name] + end var_val = get_system_state_value(state, ActivePowerVariable(), comp_type)[comp_name] if !isapprox(status_val, 0.0; atol = ABSOLUTE_TOLERANCE) - comp = get_component(ic) min = PSY.get_active_power_limits(comp).min max = PSY.get_active_power_limits(comp).max if var_val <= max && var_val >= min @@ -160,6 +165,7 @@ function update_initial_conditions!( }, } for ic in ics + isnothing(get_value(ic)) && continue var_val = get_system_state_value(state, OnVariable(), get_component_type(ic)) set_ic_quantity!(ic, var_val[get_component_name(ic)]) end From 98a66a7ae84aca5f88f957ea9b2630671deabd28 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 09:18:12 -0700 Subject: [PATCH 61/68] update obj function --- test/test_network_constructors.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 28e40011c8..2b90bb6480 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -204,7 +204,7 @@ end test_obj_values = IdDict{System, Float64}( c_sys5 => 342000.0, c_sys14 => 142000.0, - c_sys14_dc => 143000.0, + c_sys14_dc => 135000.0, ) for (ix, sys) in enumerate(systems) template = get_thermal_dispatch_template_network(DCPPowerModel) From 9bec27746af208c49f42e377cf868287d015e1f5 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 12:15:36 -0700 Subject: [PATCH 62/68] update multiplier HVDC --- src/devices_models/devices/TwoTerminalDC_branches.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/devices_models/devices/TwoTerminalDC_branches.jl b/src/devices_models/devices/TwoTerminalDC_branches.jl index 53aee2bae9..29685dbb8a 100644 --- a/src/devices_models/devices/TwoTerminalDC_branches.jl +++ b/src/devices_models/devices/TwoTerminalDC_branches.jl @@ -35,7 +35,7 @@ get_variable_multiplier( ::FlowActivePowerToFromVariable, ::Type{<:PSY.TwoTerminalHVDCLine}, ::AbstractTwoTerminalDCLineFormulation, -) = 1.0 +) = -1.0 function get_variable_multiplier( ::HVDCLosses, From 7721814cfa8e9df58882ae3b9a54f62011cb0b91 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 12:15:49 -0700 Subject: [PATCH 63/68] update tests --- test/test_device_branch_constructors.jl | 10 +++++----- test/test_network_constructors.jl | 2 +- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/test/test_device_branch_constructors.jl b/test/test_device_branch_constructors.jl index 2dc350e14f..289b2915a9 100644 --- a/test/test_device_branch_constructors.jl +++ b/test/test_device_branch_constructors.jl @@ -312,12 +312,12 @@ end available = true, active_power_flow = 0.0, # Force the flow in the opposite direction for testing purposes - active_power_limits_from = (min = -2.0, max = -2.0), - active_power_limits_to = (min = -3.0, max = 2.0), + active_power_limits_from = (min = -2.0, max = 2.0), + active_power_limits_to = (min = -2.0, max = 2.0), reactive_power_limits_from = (min = -1.0, max = 1.0), reactive_power_limits_to = (min = -1.0, max = 1.0), arc = get_arc(line), - loss = (l0 = 0.00, l1 = 0.00), + loss = LinearCurve(0.0), ) add_component!(sys_5, hvdc) @@ -331,7 +331,7 @@ end set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment) set_device_model!(template_uc, RenewableDispatch, FixedOutput) set_device_model!(template_uc, PowerLoad, StaticPowerLoad) - set_device_model!(template_uc, DeviceModel(Line, StaticBranchUnbounded)) + set_device_model!(template_uc, DeviceModel(Line, StaticBranchBounds)) set_device_model!( template_uc, DeviceModel(TwoTerminalHVDCLine, HVDCTwoTerminalLossless), @@ -424,7 +424,7 @@ end @test all( isapprox.( hvdc_ft_no_loss_values[!, "1"], - hvdc_tf_no_loss_values[!, "1"]; + -hvdc_tf_no_loss_values[!, "1"]; atol = 1e-3, ), ) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 2b90bb6480..28e40011c8 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -204,7 +204,7 @@ end test_obj_values = IdDict{System, Float64}( c_sys5 => 342000.0, c_sys14 => 142000.0, - c_sys14_dc => 135000.0, + c_sys14_dc => 143000.0, ) for (ix, sys) in enumerate(systems) template = get_thermal_dispatch_template_network(DCPPowerModel) From a8899e32a62f83cb468ecf999bf26149430aaf01 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 12:30:14 -0700 Subject: [PATCH 64/68] update test values --- test/test_network_constructors.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_network_constructors.jl b/test/test_network_constructors.jl index 28e40011c8..2b90bb6480 100644 --- a/test/test_network_constructors.jl +++ b/test/test_network_constructors.jl @@ -204,7 +204,7 @@ end test_obj_values = IdDict{System, Float64}( c_sys5 => 342000.0, c_sys14 => 142000.0, - c_sys14_dc => 143000.0, + c_sys14_dc => 135000.0, ) for (ix, sys) in enumerate(systems) template = get_thermal_dispatch_template_network(DCPPowerModel) From 5848f7755c579136b0d6361f259edd8e9fadcc74 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 16:45:01 -0700 Subject: [PATCH 65/68] add vom cost to objective --- .../common/objective_function/common.jl | 34 +++++++++++++++++++ 1 file changed, 34 insertions(+) diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index c2da8f5a3a..750556f51e 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -11,6 +11,7 @@ function add_variable_cost!( for d in devices op_cost_data = PSY.get_operation_cost(d) _add_variable_cost_to_objective!(container, U(), d, op_cost_data, V()) + _add_vom_cost_to_objective!(container, U(), d, op_cost_data, V()) end return end @@ -60,6 +61,39 @@ function add_proportional_cost!( return end +################################## +########## VOM Cost ############## +################################## + +function _add_vom_cost_to_objective!( + container::OptimizationContainer, + ::T, + component::PSY.Component, + op_cost::PSY.OperationalCost, + ::AbstractDeviceFormulation, +) where {T <: VariableType} + variable_cost = PSY.get_variable(op_cost) + power_units = PSY.get_power_units(variable_cost) + vom_cost = PSY.get_vom_cost(variable_cost) + multiplier = 1.0 # VOM Cost is always positive + cost_term = PSY.get_proportional_term(vom_cost) + iszero(cost_term) && return + base_power = get_base_power(container) + device_base_power = PSY.get_base_power(component) + cost_term_normalized = get_proportional_cost_per_system_unit( + cost_term, + power_units, + base_power, + device_base_power, + ) + for t in get_time_steps(container) + exp = + _add_proportional_term!(container, T(), d, cost_term_normalized * multiplier, t) + add_to_expression!(container, ProductionCostExpression, exp, d, t) + end + return +end + ################################## ######## OnVariable Cost ######### ################################## From 94a4e75080b4ad77dec5420a5d26827c20b578a1 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Mon, 7 Oct 2024 16:45:15 -0700 Subject: [PATCH 66/68] update name to onvar cost --- .../devices/thermal_generation.jl | 32 ++++++++++--------- 1 file changed, 17 insertions(+), 15 deletions(-) diff --git a/src/devices_models/devices/thermal_generation.jl b/src/devices_models/devices/thermal_generation.jl index ae80274519..1c99149f83 100644 --- a/src/devices_models/devices/thermal_generation.jl +++ b/src/devices_models/devices/thermal_generation.jl @@ -76,11 +76,7 @@ initial_condition_variable(::InitialTimeDurationOff, d::PSY.ThermalGen, ::Abstra ########################Objective Function################################################## # TODO: Decide what is the cost for OnVariable, if fixed or constant term in variable function proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractThermalFormulation) - return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + PSY.get_fixed(cost) -end - -function proportional_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, T::PSY.ThermalGen, U::AbstractCompactUnitCommitment) - return no_load_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + PSY.get_fixed(cost) + return onvar_cost(cost, S, T, U) + PSY.get_constant_term(PSY.get_vom_cost(PSY.get_variable(cost))) + PSY.get_fixed(cost) end proportional_cost(cost::PSY.MarketBidCost, ::OnVariable, ::PSY.ThermalGen, ::AbstractThermalFormulation) = PSY.get_no_load_cost(cost) @@ -115,17 +111,16 @@ variable_cost(cost::PSY.OperationalCost, ::PowerAboveMinimumVariable, ::PSY.Ther """ Theoretical Cost at power output zero. Mathematically is the intercept with the y-axis """ -function no_load_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, d::PSY.ThermalGen, U::AbstractThermalFormulation) - return _no_load_cost(PSY.get_variable(cost), d) +function onvar_cost(cost::PSY.ThermalGenerationCost, S::OnVariable, d::PSY.ThermalGen, U::AbstractThermalFormulation) + return _onvar_cost(PSY.get_variable(cost), d) end -function _no_load_cost(cost_function::PSY.CostCurve{PSY.PiecewisePointCurve}, d::PSY.ThermalGen) - value_curve = PSY.get_value_curve(cost_function) - cost = PSY.get_function_data(value_curve) - return last(first(PSY.get_points(cost))) +function _onvar_cost(cost_function::PSY.CostCurve{PSY.PiecewisePointCurve}, d::PSY.ThermalGen) + # OnVariableCost is included in the Point itself for PiecewisePointCurve + return 0.0 end -function _no_load_cost(cost_function::Union{PSY.CostCurve{PSY.LinearCurve}, PSY.CostCurve{PSY.QuadraticCurve}}, d::PSY.ThermalGen) +function _onvar_cost(cost_function::Union{PSY.CostCurve{PSY.LinearCurve}, PSY.CostCurve{PSY.QuadraticCurve}}, d::PSY.ThermalGen) value_curve = PSY.get_value_curve(cost_function) cost_component = PSY.get_function_data(value_curve) # Always in \$/h @@ -133,15 +128,22 @@ function _no_load_cost(cost_function::Union{PSY.CostCurve{PSY.LinearCurve}, PSY. return constant_term end -function _no_load_cost(cost_function::PSY.FuelCurve{PSY.PiecewisePointCurve}, d::PSY.ThermalGen) +function _onvar_cost(cost_function::PSY.CostCurve{PSY.PiecewiseIncrementalCurve}, d::PSY.ThermalGen) + # Input at min is used to transform to InputOutputCurve + return 0.0 +end + +function _onvar_cost(cost_function::PSY.FuelCurve{PSY.PiecewisePointCurve}, d::PSY.ThermalGen) + # OnVariableCost is included in the Point itself for PiecewisePointCurve return 0.0 end -function _no_load_cost(cost_function::PSY.FuelCurve{PSY.PiecewiseIncrementalCurve}, d::PSY.ThermalGen) +function _onvar_cost(cost_function::PSY.FuelCurve{PSY.PiecewiseIncrementalCurve}, d::PSY.ThermalGen) + # Input at min is used to transform to InputOutputCurve return 0.0 end -function _no_load_cost(cost_function::Union{PSY.FuelCurve{PSY.LinearCurve}, PSY.FuelCurve{PSY.QuadraticCurve}}, d::PSY.ThermalGen) +function _onvar_cost(cost_function::Union{PSY.FuelCurve{PSY.LinearCurve}, PSY.FuelCurve{PSY.QuadraticCurve}}, d::PSY.ThermalGen) value_curve = PSY.get_value_curve(cost_function) cost_component = PSY.get_function_data(value_curve) # In Unit/h (unit typically in ) From 7a21da9e33da2c0f45b9a523bc903fcc7b32eb00 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 8 Oct 2024 09:42:33 -0700 Subject: [PATCH 67/68] use variable cost method --- .../devices/common/objective_function/common.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/devices_models/devices/common/objective_function/common.jl b/src/devices_models/devices/common/objective_function/common.jl index 750556f51e..c02b55a51f 100644 --- a/src/devices_models/devices/common/objective_function/common.jl +++ b/src/devices_models/devices/common/objective_function/common.jl @@ -70,11 +70,11 @@ function _add_vom_cost_to_objective!( ::T, component::PSY.Component, op_cost::PSY.OperationalCost, - ::AbstractDeviceFormulation, -) where {T <: VariableType} - variable_cost = PSY.get_variable(op_cost) - power_units = PSY.get_power_units(variable_cost) - vom_cost = PSY.get_vom_cost(variable_cost) + ::U, +) where {T <: VariableType, U <: AbstractDeviceFormulation} + variable_cost_data = variable_cost(op_cost, T(), component, U()) + power_units = PSY.get_power_units(variable_cost_data) + vom_cost = PSY.get_vom_cost(variable_cost_data) multiplier = 1.0 # VOM Cost is always positive cost_term = PSY.get_proportional_term(vom_cost) iszero(cost_term) && return From 1482504b7d72c711d0b5c8148b1383e6d5e471b6 Mon Sep 17 00:00:00 2001 From: rodrigomha Date: Tue, 8 Oct 2024 09:52:39 -0700 Subject: [PATCH 68/68] remove old constraint export --- docs/src/api/PowerSimulations.md | 1 - 1 file changed, 1 deletion(-) diff --git a/docs/src/api/PowerSimulations.md b/docs/src/api/PowerSimulations.md index be88337f83..a00e3bce49 100644 --- a/docs/src/api/PowerSimulations.md +++ b/docs/src/api/PowerSimulations.md @@ -269,7 +269,6 @@ FlowLimitConstraint FlowRateConstraint FlowRateConstraintFromTo FlowRateConstraintToFrom -HVDCLossesAbsoluteValue HVDCPowerBalance NetworkFlowConstraint RateLimitConstraint