Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jd/implement state estimator #5

Open
wants to merge 33 commits into
base: main
Choose a base branch
from
Open

Conversation

jd-lara
Copy link
Member

@jd-lara jd-lara commented Jul 10, 2024

Thanks for opening a PR to PowerSimulationsDecomposition.jl, please take note of the following when making a PR:

Check the contributor guidelines

@jd-lara jd-lara requested a review from rodrigomha July 10, 2024 21:45
@jd-lara jd-lara self-assigned this Jul 10, 2024
Copy link

@github-actions github-actions bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

selected_line=["CA-1", "CB-1", "AB1"]


[JuliaFormatter] reported by reviewdog 🐶

monitoredlined_line=["A28"]
limit=0.2


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,se=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()


[JuliaFormatter] reported by reviewdog 🐶

bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[], PowerLoad=[], ActivePowerBalance__ACBus=[],
StateEstimationInjections__ACBus=[], ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn], bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], ptdf[line,bn],
bus_ActivePowerBalance__ACBus[t,bn]*ptdf[line,bn],bus_StateEstimationInjections__ACBus[t,bn]*ptdf[line,bn]))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function check_models(uc0,uc2a)
sol0=Dict()
v0=uc0.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol0=Dict(zip(name.(v0),value.(v0)))


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end


[JuliaFormatter] reported by reviewdog 🐶

v01=uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")]


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end
v2a=uc2a.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol2a=Dict(zip(name.(v2a),value.(v2a)))
for (k,v) in sol2a
if abs(sol2a[k]-sol0[k])>0.00001
println("k,",k,",",sol2a[k],",",sol0[k])
var0=variable_by_name(uc0.JuMPmodel,k)
fix(var0,sol2a[k],force=true)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶

end
return(objd0,sol0, objv0)
end


[JuliaFormatter] reported by reviewdog 🐶

name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "2"),


[JuliaFormatter] reported by reviewdog 🐶

name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "2"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "2"),


[JuliaFormatter] reported by reviewdog 🐶

name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "3"),


[JuliaFormatter] reported by reviewdog 🐶

name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "2"), to_area=get_component(Area, sys2, "3"),


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc, DeviceModel(MonitoredLine, StaticBranchUnbounded, use_slacks = true))


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template_uc2,
DeviceModel(Line, StaticBranchUnbounded; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template_uc2, DeviceModel(MonitoredLine, StaticBranch, use_slacks = true))


[JuliaFormatter] reported by reviewdog 🐶

set_flow_limits!(l,(from_to=limit,to_from=limit))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)


[JuliaFormatter] reported by reviewdog 🐶

feedforwards=Dict(
"UC_Subsystem" => uc_simulation_ff,
),


[JuliaFormatter] reported by reviewdog 🐶

uc0=models.decision_models[1].internal.container
uc2=models.decision_models[2].internal.container
println("obj uc0,uc2,",objective_value(uc0.JuMPmodel),",",objective_value(uc2.JuMPmodel),",",objective_value(uc2b.JuMPmodel),",diff,",objective_value(uc0.JuMPmodel)-objective_value(uc2.JuMPmodel))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for k in all_constraints(uc0.JuMPmodel,; include_variable_in_set_constraints = true)
println(name(k),",",k)
end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for k in all_constraints(uc2.JuMPmodel,; include_variable_in_set_constraints = true)
println(name(k),",",k)
end


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in uc0.constraints println(k) end
uc0.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,1]
for (k,v) in uc2.constraints println(k) end
uc2.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,1]
value.(uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")][:,1])
fix_value.(uc2.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")][:,1])


[JuliaFormatter] reported by reviewdog 🐶

bus_df_uc0=read_bus_df(results_uc0,0)
bus_df_uc2=read_bus_df(results_uc1,0)
bus_df_uc0[!,"buscheck"]=bus_df_uc0[!,"ActivePowerBalance__ACBus"]-bus_df_uc0[!,"ThermalStandard"]/100-bus_df_uc0[!,"PowerLoad"]/100
println("UC0 bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df_uc0))
bus_df_uc2[!,"buscheck"]=bus_df_uc2[!,"ActivePowerBalance__ACBus"]-bus_df_uc2[!,"ThermalStandard"]/100-bus_df_uc2[!,"PowerLoad"]/100
println("UC2 bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df_uc2))


[JuliaFormatter] reported by reviewdog 🐶

area_df0=combine(groupby(bus_df_uc0, [:t, :area]), [:ThermalStandard, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum)
area_df0[!,"areacheck"]=area_df0[!,"ActivePowerBalance__ACBus_sum"]-area_df0[!,"ThermalStandard_sum"]/100-area_df0[!,"PowerLoad_sum"]/100
sum(area_df0[!,"areacheck"])
area1_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), area_df0)
area2_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), area_df0)
area3_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), area_df0)
area_df2=combine(groupby(bus_df_uc2, [:t, :area]), [:ThermalStandard, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum)
area_df2[!,"areacheck"]=area_df2[!,"ActivePowerBalance__ACBus_sum"]-area_df2[!,"ThermalStandard_sum"]/100-area_df2[!,"PowerLoad_sum"]/100
sum(area_df2[!,"areacheck"])
area1_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), area_df2)
area2_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), area_df2)
area3_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), area_df2)
println("area1 check uc0-uc2,", sum(area1_df2[!,"ActivePowerBalance__ACBus_sum"]-area1_df0[!,"ActivePowerBalance__ACBus_sum"]))
println("area2 check uc0-uc2,", sum(area2_df2[!,"ActivePowerBalance__ACBus_sum"]-area2_df0[!,"ActivePowerBalance__ACBus_sum"]))
println("area3 check uc0-uc2,", sum(area3_df2[!,"ActivePowerBalance__ACBus_sum"]-area3_df0[!,"ActivePowerBalance__ACBus_sum"]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow_df0=combine(groupby(bus_df_uc0, [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum)
areaflow1_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), areaflow_df0)
areaflow2_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), areaflow_df0)
areaflow3_df0=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), areaflow_df0)
areaflow_df2=combine(groupby(bus_df_uc2, [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum)
areaflow1_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), areaflow_df2)
areaflow2_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), areaflow_df2)
areaflow3_df2=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), areaflow_df2)


[JuliaFormatter] reported by reviewdog 🐶

println("areaflow1 check sum(uc0 flowcontribution_sum - uc2 flowcontribution_sum,",sum(areaflow1_df0[!,"flowcontribution_sum"]-areaflow1_df2[!,"flowcontribution_sum"]))
println("areaflow2 check sum(uc0 flowcontribution_sum - uc2 flowcontribution_sum,",sum(areaflow2_df0[!,"flowcontribution_sum"]-areaflow2_df2[!,"flowcontribution_sum"]))
println("areaflow3 check sum(uc0 flowcontribution_sum - uc2 flowcontribution_sum,",sum(areaflow3_df0[!,"flowcontribution_sum"]-areaflow3_df2[!,"flowcontribution_sum"]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow1_df0[!,"seflow"]=areaflow1_df0[!,"flowcontribution_sum"]+areaflow2_df0[!,"flowcontribution_sum"]+areaflow3_df0[!,"flowcontribution_sum"]
areaflow1_df2[!,"seflow"]=areaflow1_df2[!,"flowcontribution_sum"]+areaflow2_df2[!,"flowcontribution_sum"]+areaflow3_df2[!,"flowcontribution_sum"]


[JuliaFormatter] reported by reviewdog 🐶

println("varflow_uc0,",value.(uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, MonitoredLine}("")]["A28",:]) )
println("varflow_uc2,",value.(uc2.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, MonitoredLine}("")]["A28",:]) )


[JuliaFormatter] reported by reviewdog 🐶

NT=24
selected_line=["CA-1", "CB-1", "AB1"]
monitoredlined_line=["A28"]
limit=0.2
use_monitoredline=1
sys = build_system(PSISystems, "modified_RTS_GMLC_DA_sys")
sys2 = build_system(PSISystems, "modified_RTS_GMLC_DA_sys")
transform_single_time_series!(sys, Hour(NT), Hour(NT))
transform_single_time_series!(sys2, Hour(1), Hour(1))
#transform_single_time_series!(sys, Hour(1), Hour(1))
#transform_single_time_series!(sys2, Hour(1), Hour(1))
#enforced_region = ["1", "2", "3"] # first one is MRTO, second one is M2M coordinator, third one in RT is loopflow under M2M
enforced_region=["1","3","2"]
#enforced_region=["2","3","1"]
ptdf=PTDF(sys)
buildsubsystem(sys2, enforced_region, union(selected_line,monitoredlined_line))
exchange_1_2 = AreaInterchange(;
name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "2"),
flow_limits=(from_to=99999, to_from=99999),
)
add_component!(sys, exchange_1_2)
exchange_1_3 = AreaInterchange(;
name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "3"),
flow_limits=(from_to=99999, to_from=99999),
)
add_component!(sys, exchange_1_3)
exchange_2_3 = AreaInterchange(;
name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "2"), to_area=get_component(Area, sys, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys, exchange_2_3)
exchange_1_2 = AreaInterchange(;
name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "2"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_2)
add_component_to_subsystem!(sys2, "a", exchange_1_2)
add_component_to_subsystem!(sys2, "b", exchange_1_2)


[JuliaFormatter] reported by reviewdog 🐶

exchange_1_3 = AreaInterchange(;
name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_3)
add_component_to_subsystem!(sys2, "a", exchange_1_3)
add_component_to_subsystem!(sys2, "b", exchange_1_3)
exchange_2_3 = AreaInterchange(;
name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "2"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_2_3)
add_component_to_subsystem!(sys2, "a", exchange_2_3)
add_component_to_subsystem!(sys2, "b", exchange_2_3)
template_uc = ProblemTemplate(NetworkModel(AreaPTDFPowerModel; use_slacks=true))
set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc, PowerLoad, StaticPowerLoad)
set_device_model!(template_uc, AreaInterchange, StaticBranch)
if use_monitoredline==0
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),
#DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
#DeviceModel(Line, StaticBranch; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
))
else
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys, b)
PSY.convert_component!(sys, line, MonitoredLine)
end
end


[JuliaFormatter] reported by reviewdog 🐶

template_uc2 = MultiProblemTemplate(NetworkModel(SplitAreaPTDFPowerModel; use_slacks=true), ["a", "b"])
set_device_model!(template_uc2, AreaInterchange, StaticBranch)
set_device_model!(template_uc2, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc2, PowerLoad, StaticPowerLoad)
if use_monitoredline==0
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranch; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
l = PSY.get_component(ACBranch, sys2, "A28")
set_rating!(l,limit)
else
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranchUnbounded; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc2, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys2, b)
PSY.convert_component!(sys2, line, MonitoredLine)
end
l = PSY.get_component(MonitoredLine, sys2, "A28")
set_flow_limits!(l,(from_to=limit,to_from=limit))


[JuliaFormatter] reported by reviewdog 🐶

models = SimulationModels(;
decision_models=[
DecisionModel(
template_uc,
sys;
name="UC0",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,


[JuliaFormatter] reported by reviewdog 🐶

),
system_to_file=false,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

DecisionModel(
MultiRegionProblem,
template_uc2,
sys2;
name="UC_Subsystem",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,
#Xpress.Optimizer,
#"MIPRELSTOP" => 0.00, # Set the relative mip gap tolerance
#"MAXMEMORYSOFT" => 600000, # Set the maximum amount of memory the solver can use (in MB)
),
system_to_file=false,
initialize_model=true,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
rebuild_model=false,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

],
)
uc_simulation_ff = Vector{PowerSimulations.AbstractAffectFeedforward}()
#FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)
#push!(uc_simulation_ff, FVFF_area_interchange)


[JuliaFormatter] reported by reviewdog 🐶

#LBFF = FixValueFeedforward(;
# component_type=ThermalStandard,
# source=OnVariable,
# affected_values=[OnVariable],
#)
#push!(uc_simulation_ff, LBFF)


[JuliaFormatter] reported by reviewdog 🐶

sequence = SimulationSequence(;
models=models,
feedforwards=Dict(
"UC_Subsystem" => uc_simulation_ff,
),
ini_cond_chronology=InterProblemChronology(),
);
# use different names for saving the solution
sim = Simulation(;
name="sim",
steps=1,
models=models,
sequence=sequence,
initial_time=DateTime("2020-01-01T00:00:00"),
simulation_folder=mktempdir(),
);
build_out = build!(sim; console_level=Logging.Info, serialize=false)
execute_status = execute!(sim; enable_progress_bar=true);
uc0=models.decision_models[1].internal.container
uc2b=models.decision_models[2].internal.container.subproblems["b"]
uc2a=models.decision_models[2].internal.container.subproblems["a"]
println("obj uc0,uc2a,uc2b,",objective_value(uc0.JuMPmodel),",",objective_value(uc2a.JuMPmodel),",",objective_value(uc2b.JuMPmodel),",diff,",objective_value(uc0.JuMPmodel)-objective_value(uc2a.JuMPmodel)-objective_value(uc2b.JuMPmodel))
results = SimulationResults(sim)
# #results = SimulationResults(sim)
results_uc = get_decision_problem_results(results, "UC0")
results_ha = get_decision_problem_results(results, "UC_Subsystem")
aa = read_realized_variable(results_uc, "FlowActivePowerVariable__MonitoredLine")
ac = read_realized_variable(results_ha, "FlowActivePowerVariable__MonitoredLine")
area_df=DataFrames.DataFrame(;
t =Int64[], uc=[], area1_gen = [], area1_load = [], area2_gen = [], area2_load = [],
area3_gen = [], area3_load = [],


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

area_df=combine(groupby(bus_df, [:t, :area]), [:ThermalStandard, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum)
area_df[!,"areacheck"]=area_df[!,"ActivePowerBalance__ACBus_sum"]-area_df[!,"ThermalStandard_sum"]/100-area_df[!,"PowerLoad_sum"]/100
area1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), area_df)
area2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), area_df)
area3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), area_df)
area1_df[!,"secheck"]=lag(area1_df[!,"ActivePowerBalance__ACBus_sum"])-area1_df[!,"StateEstimationInjections__ACBus_sum"]
area2_df[!,"secheck"]=lag(area2_df[!,"ActivePowerBalance__ACBus_sum"])-area2_df[!,"StateEstimationInjections__ACBus_sum"]
area3_df[!,"secheck"]=lag(area3_df[!,"ActivePowerBalance__ACBus_sum"])-area3_df[!,"StateEstimationInjections__ACBus_sum"]
println("area1 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("area2 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("area3 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow_df=combine(groupby(bus_df, [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum)


[JuliaFormatter] reported by reviewdog 🐶

areaflow1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), areaflow_df)
areaflow2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), areaflow_df)
areaflow3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), areaflow_df)
areaflow1_df[!,"secheck"]=lag(areaflow1_df[!,"flowcontribution_sum"])-areaflow1_df[!,"loopflowcontribution_sum"]
areaflow2_df[!,"secheck"]=lag(areaflow2_df[!,"flowcontribution_sum"])-areaflow2_df[!,"loopflowcontribution_sum"]
areaflow3_df[!,"secheck"]=lag(areaflow3_df[!,"flowcontribution_sum"])-areaflow3_df[!,"loopflowcontribution_sum"]
println("areaflow1 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("areaflow2 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("areaflow3 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow13_df=innerjoin(areaflow1_df, areaflow3_df, on=:t,renamecols = "_1" => "_3")
areaflow132_df=innerjoin(areaflow13_df, areaflow2_df, on=:t,renamecols = "" => "_2")
areaflow132_df[!,"flow_subsystem_a"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"loopflowcontribution_sum_2"]
areaflow132_df[!,"flow_subsystem_b"]=areaflow132_df[!,"loopflowcontribution_sum_1"]+areaflow132_df[!,"loopflowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
areaflow132_df[!,"seflow"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
StatsPlots.@df(areaflow132_df,Plots.plot(:t,[:flow_subsystem_a :flow_subsystem_b :seflow]))


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,line,se=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()


[JuliaFormatter] reported by reviewdog 🐶

bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[], PowerLoad=[], ActivePowerBalance__ACBus=[],
StateEstimationInjections__ACBus=[], ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn], bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], ptdf[line,bn],
bus_ActivePowerBalance__ACBus[t,bn]*ptdf[line,bn],bus_StateEstimationInjections__ACBus[t,bn]*ptdf[line,bn]))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function check_models(uc0,uc2a)
sol0=Dict()
v0=uc0.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol0=Dict(zip(name.(v0),value.(v0)))


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end


[JuliaFormatter] reported by reviewdog 🐶

v01=uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")]


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end
v2a=uc2a.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol2a=Dict(zip(name.(v2a),value.(v2a)))


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in sol2a
if abs(sol2a[k]-sol0[k])>0.00001
println("k,",k,",",sol2a[k],",",sol0[k])
var0=variable_by_name(uc0.JuMPmodel,k)
fix(var0,sol2a[k],force=true)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶

end
return(objd0,sol0, objv0)
end


[JuliaFormatter] reported by reviewdog 🐶

function write_param_acbus(name,uc2a,line)
parameter_ACBus_2a=uc2a.parameters[InfrastructureSystems.Optimization.ParameterKey{PowerSimulationsDecomposition.StateEstimationInjections, ACBus}("")]
i=axes(parameter_ACBus_2a.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_ACBus_2a.parameter_array[j1,j2]
println(name,",j1,j2,",j1,",",j2,",value,",value(parameter_ACBus_2a.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_coeff(uc2a,linename,t,fnm,use_monitoredline)
fixedmw=0; fixedflow=0; areamw=Dict(); areaflow=Dict(); areaload=Dict(); arealoadflow=Dict()
for j=1:3 areamw[j]=0; areaflow[j]=0 end
c=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,t]


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(c[i], k)!=0


[JuliaFormatter] reported by reviewdog 🐶

if name(k)=="param"
areaload[i]=-normalized_coefficient(c[i], k)
println("area,",i,",",areaload[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if use_monitoredline==1
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][linename,t]
else
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, Line}("")][linename,t]
end
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
if name(k)=="param"
fixedmw=fixedmw+value(k);
fixedflow=fixedflow+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶

str=split(name(k),"{");
println("name,",name(k),",str,",str)
if first(str[1],1)=="A"
a= parse(Int64,first(str[2],1));# tt=firstindex(split(str[2],", "))
areamw[a]=areamw[a]+value(k)
areaflow[a]=areaflow[a]+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

return(fixedmw,areamw,fixedflow,areaflow)


[JuliaFormatter] reported by reviewdog 🐶

function write_param_internalload(name,uc0,line)
parameter_load_0=uc0.parameters[InfrastructureSystems.Optimization.ParameterKey{ActivePowerTimeSeriesParameter, PowerLoad}("")]
i=axes(parameter_load_0.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc0.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_load_0.parameter_array[j1,j2]
println(name,",load,j1,j2,",j1,",",j2,",value,",value(parameter_load_0.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end
function write_acbus(ucname,uc0,ax,t, debug=0)
area_param=Dict(); area_gen=Dict()
area_param[1]=0;area_param[2]=0;area_param[3]=0
area_gen[1]=0;area_gen[2]=0;area_gen[3]=0
for bus in ax[1]
for k in keys(PSI.get_expression(uc0, ActivePowerBalance(), PSY.ACBus)[bus,t].terms)
coef=PSI.get_expression(uc0, ActivePowerBalance(),PSY.ACBus)[bus,t].terms[k]
v=value(k);
if debug==1
println(ucname,";bus;t;",bus,";",t,";",k,";coefficient;", coef,";solution;",v,";coef*v;",coef*v)
end
n=Int64(floor(bus/100))
if name(k)=="param" area_param[n]=area_param[n]+coef*v else area_gen[n]=area_gen[n]+coef*v end


[JuliaFormatter] reported by reviewdog 🐶

end
return(area_gen,area_param)


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(get_bus(b))) == enforced_region[1] || get_name(get_area(get_bus(b))) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(b)) == enforced_region[1] || get_name(get_area(b)) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, "b", l)
end
end

Copy link

@github-actions github-actions bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys, exchange_2_3)
exchange_1_2 = AreaInterchange(;
name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "2"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_2)
add_component_to_subsystem!(sys2, "a", exchange_1_2)
add_component_to_subsystem!(sys2, "b", exchange_1_2)


[JuliaFormatter] reported by reviewdog 🐶

exchange_1_3 = AreaInterchange(;
name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_3)
add_component_to_subsystem!(sys2, "a", exchange_1_3)
add_component_to_subsystem!(sys2, "b", exchange_1_3)
exchange_2_3 = AreaInterchange(;
name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "2"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_2_3)
add_component_to_subsystem!(sys2, "a", exchange_2_3)
add_component_to_subsystem!(sys2, "b", exchange_2_3)
template_uc = ProblemTemplate(NetworkModel(AreaPTDFPowerModel; use_slacks=true))
set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc, PowerLoad, StaticPowerLoad)
set_device_model!(template_uc, AreaInterchange, StaticBranch)
if use_monitoredline==0
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),
#DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
#DeviceModel(Line, StaticBranch; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
))
else
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys, b)
PSY.convert_component!(sys, line, MonitoredLine)
end
end


[JuliaFormatter] reported by reviewdog 🐶

template_uc2 = MultiProblemTemplate(NetworkModel(SplitAreaPTDFPowerModel; use_slacks=true), ["a", "b"])
set_device_model!(template_uc2, AreaInterchange, StaticBranch)
set_device_model!(template_uc2, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc2, PowerLoad, StaticPowerLoad)
if use_monitoredline==0
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranch; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
l = PSY.get_component(ACBranch, sys2, "A28")
set_rating!(l,limit)
else
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranchUnbounded; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc2, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys2, b)
PSY.convert_component!(sys2, line, MonitoredLine)
end
l = PSY.get_component(MonitoredLine, sys2, "A28")
set_flow_limits!(l,(from_to=limit,to_from=limit))


[JuliaFormatter] reported by reviewdog 🐶

models = SimulationModels(;
decision_models=[
DecisionModel(
template_uc,
sys;
name="UC0",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,


[JuliaFormatter] reported by reviewdog 🐶

),
system_to_file=false,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

DecisionModel(
MultiRegionProblem,
template_uc2,
sys2;
name="UC_Subsystem",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,
#Xpress.Optimizer,
#"MIPRELSTOP" => 0.00, # Set the relative mip gap tolerance
#"MAXMEMORYSOFT" => 600000, # Set the maximum amount of memory the solver can use (in MB)
),
system_to_file=false,
initialize_model=true,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
rebuild_model=false,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

],
)
uc_simulation_ff = Vector{PowerSimulations.AbstractAffectFeedforward}()
#FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)
#push!(uc_simulation_ff, FVFF_area_interchange)


[JuliaFormatter] reported by reviewdog 🐶

#LBFF = FixValueFeedforward(;
# component_type=ThermalStandard,
# source=OnVariable,
# affected_values=[OnVariable],
#)
#push!(uc_simulation_ff, LBFF)


[JuliaFormatter] reported by reviewdog 🐶

sequence = SimulationSequence(;
models=models,
feedforwards=Dict(
"UC_Subsystem" => uc_simulation_ff,
),
ini_cond_chronology=InterProblemChronology(),
);
# use different names for saving the solution
sim = Simulation(;
name="sim",
steps=1,
models=models,
sequence=sequence,
initial_time=DateTime("2020-01-01T00:00:00"),
simulation_folder=mktempdir(),
);
build_out = build!(sim; console_level=Logging.Info, serialize=false)
execute_status = execute!(sim; enable_progress_bar=true);
uc0=models.decision_models[1].internal.container
uc2b=models.decision_models[2].internal.container.subproblems["b"]
uc2a=models.decision_models[2].internal.container.subproblems["a"]
println("obj uc0,uc2a,uc2b,",objective_value(uc0.JuMPmodel),",",objective_value(uc2a.JuMPmodel),",",objective_value(uc2b.JuMPmodel),",diff,",objective_value(uc0.JuMPmodel)-objective_value(uc2a.JuMPmodel)-objective_value(uc2b.JuMPmodel))
results = SimulationResults(sim)
# #results = SimulationResults(sim)
results_uc = get_decision_problem_results(results, "UC0")
results_ha = get_decision_problem_results(results, "UC_Subsystem")
aa = read_realized_variable(results_uc, "FlowActivePowerVariable__MonitoredLine")
ac = read_realized_variable(results_ha, "FlowActivePowerVariable__MonitoredLine")
area_df=DataFrames.DataFrame(;
t =Int64[], uc=[], area1_gen = [], area1_load = [], area2_gen = [], area2_load = [],
area3_gen = [], area3_load = [],


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

area_df=combine(groupby(bus_df, [:t, :area]), [:ThermalStandard, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum)
area_df[!,"areacheck"]=area_df[!,"ActivePowerBalance__ACBus_sum"]-area_df[!,"ThermalStandard_sum"]/100-area_df[!,"PowerLoad_sum"]/100
area1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), area_df)
area2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), area_df)
area3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), area_df)
area1_df[!,"secheck"]=lag(area1_df[!,"ActivePowerBalance__ACBus_sum"])-area1_df[!,"StateEstimationInjections__ACBus_sum"]
area2_df[!,"secheck"]=lag(area2_df[!,"ActivePowerBalance__ACBus_sum"])-area2_df[!,"StateEstimationInjections__ACBus_sum"]
area3_df[!,"secheck"]=lag(area3_df[!,"ActivePowerBalance__ACBus_sum"])-area3_df[!,"StateEstimationInjections__ACBus_sum"]
println("area1 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("area2 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("area3 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow_df=combine(groupby(bus_df, [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum)


[JuliaFormatter] reported by reviewdog 🐶

areaflow1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), areaflow_df)
areaflow2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), areaflow_df)
areaflow3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), areaflow_df)
areaflow1_df[!,"secheck"]=lag(areaflow1_df[!,"flowcontribution_sum"])-areaflow1_df[!,"loopflowcontribution_sum"]
areaflow2_df[!,"secheck"]=lag(areaflow2_df[!,"flowcontribution_sum"])-areaflow2_df[!,"loopflowcontribution_sum"]
areaflow3_df[!,"secheck"]=lag(areaflow3_df[!,"flowcontribution_sum"])-areaflow3_df[!,"loopflowcontribution_sum"]
println("areaflow1 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("areaflow2 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("areaflow3 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow13_df=innerjoin(areaflow1_df, areaflow3_df, on=:t,renamecols = "_1" => "_3")
areaflow132_df=innerjoin(areaflow13_df, areaflow2_df, on=:t,renamecols = "" => "_2")
areaflow132_df[!,"flow_subsystem_a"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"loopflowcontribution_sum_2"]
areaflow132_df[!,"flow_subsystem_b"]=areaflow132_df[!,"loopflowcontribution_sum_1"]+areaflow132_df[!,"loopflowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
areaflow132_df[!,"seflow"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
StatsPlots.@df(areaflow132_df,Plots.plot(:t,[:flow_subsystem_a :flow_subsystem_b :seflow]))


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,line,se=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()


[JuliaFormatter] reported by reviewdog 🐶

bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[], PowerLoad=[], ActivePowerBalance__ACBus=[],
StateEstimationInjections__ACBus=[], ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn], bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], ptdf[line,bn],
bus_ActivePowerBalance__ACBus[t,bn]*ptdf[line,bn],bus_StateEstimationInjections__ACBus[t,bn]*ptdf[line,bn]))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function check_models(uc0,uc2a)
sol0=Dict()
v0=uc0.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol0=Dict(zip(name.(v0),value.(v0)))


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end


[JuliaFormatter] reported by reviewdog 🐶

v01=uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")]


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end
v2a=uc2a.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol2a=Dict(zip(name.(v2a),value.(v2a)))


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in sol2a
if abs(sol2a[k]-sol0[k])>0.00001
println("k,",k,",",sol2a[k],",",sol0[k])
var0=variable_by_name(uc0.JuMPmodel,k)
fix(var0,sol2a[k],force=true)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶

end
return(objd0,sol0, objv0)
end


[JuliaFormatter] reported by reviewdog 🐶

function write_param_acbus(name,uc2a,line)
parameter_ACBus_2a=uc2a.parameters[InfrastructureSystems.Optimization.ParameterKey{PowerSimulationsDecomposition.StateEstimationInjections, ACBus}("")]
i=axes(parameter_ACBus_2a.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_ACBus_2a.parameter_array[j1,j2]
println(name,",j1,j2,",j1,",",j2,",value,",value(parameter_ACBus_2a.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_coeff(uc2a,linename,t,fnm,use_monitoredline)
fixedmw=0; fixedflow=0; areamw=Dict(); areaflow=Dict(); areaload=Dict(); arealoadflow=Dict()
for j=1:3 areamw[j]=0; areaflow[j]=0 end
c=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,t]


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(c[i], k)!=0


[JuliaFormatter] reported by reviewdog 🐶

if name(k)=="param"
areaload[i]=-normalized_coefficient(c[i], k)
println("area,",i,",",areaload[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if use_monitoredline==1
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][linename,t]
else
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, Line}("")][linename,t]
end
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
if name(k)=="param"
fixedmw=fixedmw+value(k);
fixedflow=fixedflow+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶

str=split(name(k),"{");
println("name,",name(k),",str,",str)
if first(str[1],1)=="A"
a= parse(Int64,first(str[2],1));# tt=firstindex(split(str[2],", "))
areamw[a]=areamw[a]+value(k)
areaflow[a]=areaflow[a]+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

return(fixedmw,areamw,fixedflow,areaflow)


[JuliaFormatter] reported by reviewdog 🐶

function write_param_internalload(name,uc0,line)
parameter_load_0=uc0.parameters[InfrastructureSystems.Optimization.ParameterKey{ActivePowerTimeSeriesParameter, PowerLoad}("")]
i=axes(parameter_load_0.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc0.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_load_0.parameter_array[j1,j2]
println(name,",load,j1,j2,",j1,",",j2,",value,",value(parameter_load_0.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end
function write_acbus(ucname,uc0,ax,t, debug=0)
area_param=Dict(); area_gen=Dict()
area_param[1]=0;area_param[2]=0;area_param[3]=0
area_gen[1]=0;area_gen[2]=0;area_gen[3]=0
for bus in ax[1]
for k in keys(PSI.get_expression(uc0, ActivePowerBalance(), PSY.ACBus)[bus,t].terms)
coef=PSI.get_expression(uc0, ActivePowerBalance(),PSY.ACBus)[bus,t].terms[k]
v=value(k);
if debug==1
println(ucname,";bus;t;",bus,";",t,";",k,";coefficient;", coef,";solution;",v,";coef*v;",coef*v)
end
n=Int64(floor(bus/100))
if name(k)=="param" area_param[n]=area_param[n]+coef*v else area_gen[n]=area_gen[n]+coef*v end


[JuliaFormatter] reported by reviewdog 🐶

end
return(area_gen,area_param)


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(get_bus(b))) == enforced_region[1] || get_name(get_area(get_bus(b))) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(b)) == enforced_region[1] || get_name(get_area(b)) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, "b", l)
end
end

Copy link

@github-actions github-actions bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys, exchange_2_3)
exchange_1_2 = AreaInterchange(;
name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "2"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_2)
add_component_to_subsystem!(sys2, "a", exchange_1_2)
add_component_to_subsystem!(sys2, "b", exchange_1_2)


[JuliaFormatter] reported by reviewdog 🐶

exchange_1_3 = AreaInterchange(;
name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "1"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_1_3)
add_component_to_subsystem!(sys2, "a", exchange_1_3)
add_component_to_subsystem!(sys2, "b", exchange_1_3)
exchange_2_3 = AreaInterchange(;
name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys2, "2"), to_area=get_component(Area, sys2, "3"),
flow_limits=(from_to=99999, to_from=99999),


[JuliaFormatter] reported by reviewdog 🐶

add_component!(sys2, exchange_2_3)
add_component_to_subsystem!(sys2, "a", exchange_2_3)
add_component_to_subsystem!(sys2, "b", exchange_2_3)
template_uc = ProblemTemplate(NetworkModel(AreaPTDFPowerModel; use_slacks=true))
set_device_model!(template_uc, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc, PowerLoad, StaticPowerLoad)
set_device_model!(template_uc, AreaInterchange, StaticBranch)
if use_monitoredline==0
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),
#DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
#DeviceModel(Line, StaticBranch; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),
))
else
set_device_model!(template_uc,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys, b)
PSY.convert_component!(sys, line, MonitoredLine)
end
end


[JuliaFormatter] reported by reviewdog 🐶

template_uc2 = MultiProblemTemplate(NetworkModel(SplitAreaPTDFPowerModel; use_slacks=true), ["a", "b"])
set_device_model!(template_uc2, AreaInterchange, StaticBranch)
set_device_model!(template_uc2, ThermalStandard, ThermalBasicUnitCommitment)
set_device_model!(template_uc2, PowerLoad, StaticPowerLoad)
if use_monitoredline==0
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranch; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
l = PSY.get_component(ACBranch, sys2, "A28")
set_rating!(l,limit)
else
set_device_model!(template_uc2,
DeviceModel(Line, StaticBranchUnbounded; use_slacks=true, attributes=Dict("filter_function" => x -> get_name(x) in union(selected_line,monitoredlined_line)),))
set_device_model!(template_uc2, MonitoredLine, StaticBranch)
for b in monitoredlined_line
line = PSY.get_component(Line, sys2, b)
PSY.convert_component!(sys2, line, MonitoredLine)
end
l = PSY.get_component(MonitoredLine, sys2, "A28")
set_flow_limits!(l,(from_to=limit,to_from=limit))


[JuliaFormatter] reported by reviewdog 🐶

models = SimulationModels(;
decision_models=[
DecisionModel(
template_uc,
sys;
name="UC0",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,


[JuliaFormatter] reported by reviewdog 🐶

),
system_to_file=false,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

DecisionModel(
MultiRegionProblem,
template_uc2,
sys2;
name="UC_Subsystem",
optimizer=optimizer_with_attributes(
HiGHS.Optimizer,
#Xpress.Optimizer,
#"MIPRELSTOP" => 0.00, # Set the relative mip gap tolerance
#"MAXMEMORYSOFT" => 600000, # Set the maximum amount of memory the solver can use (in MB)
),
system_to_file=false,
initialize_model=true,
optimizer_solve_log_print=false,
direct_mode_optimizer=true,
rebuild_model=false,
store_variable_names=true,
calculate_conflict=true,


[JuliaFormatter] reported by reviewdog 🐶

],
)
uc_simulation_ff = Vector{PowerSimulations.AbstractAffectFeedforward}()
#FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)
#push!(uc_simulation_ff, FVFF_area_interchange)


[JuliaFormatter] reported by reviewdog 🐶

#LBFF = FixValueFeedforward(;
# component_type=ThermalStandard,
# source=OnVariable,
# affected_values=[OnVariable],
#)
#push!(uc_simulation_ff, LBFF)


[JuliaFormatter] reported by reviewdog 🐶

sequence = SimulationSequence(;
models=models,
feedforwards=Dict(
"UC_Subsystem" => uc_simulation_ff,
),
ini_cond_chronology=InterProblemChronology(),
);
# use different names for saving the solution
sim = Simulation(;
name="sim",
steps=1,
models=models,
sequence=sequence,
initial_time=DateTime("2020-01-01T00:00:00"),
simulation_folder=mktempdir(),
);
build_out = build!(sim; console_level=Logging.Info, serialize=false)
execute_status = execute!(sim; enable_progress_bar=true);
uc0=models.decision_models[1].internal.container
uc2b=models.decision_models[2].internal.container.subproblems["b"]
uc2a=models.decision_models[2].internal.container.subproblems["a"]
println("obj uc0,uc2a,uc2b,",objective_value(uc0.JuMPmodel),",",objective_value(uc2a.JuMPmodel),",",objective_value(uc2b.JuMPmodel),",diff,",objective_value(uc0.JuMPmodel)-objective_value(uc2a.JuMPmodel)-objective_value(uc2b.JuMPmodel))
results = SimulationResults(sim)
# #results = SimulationResults(sim)
results_uc = get_decision_problem_results(results, "UC0")
results_ha = get_decision_problem_results(results, "UC_Subsystem")
aa = read_realized_variable(results_uc, "FlowActivePowerVariable__MonitoredLine")
ac = read_realized_variable(results_ha, "FlowActivePowerVariable__MonitoredLine")
area_df=DataFrames.DataFrame(;
t =Int64[], uc=[], area1_gen = [], area1_load = [], area2_gen = [], area2_load = [],
area3_gen = [], area3_load = [],


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

area_df=combine(groupby(bus_df, [:t, :area]), [:ThermalStandard, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum)
area_df[!,"areacheck"]=area_df[!,"ActivePowerBalance__ACBus_sum"]-area_df[!,"ThermalStandard_sum"]/100-area_df[!,"PowerLoad_sum"]/100
area1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), area_df)
area2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), area_df)
area3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), area_df)
area1_df[!,"secheck"]=lag(area1_df[!,"ActivePowerBalance__ACBus_sum"])-area1_df[!,"StateEstimationInjections__ACBus_sum"]
area2_df[!,"secheck"]=lag(area2_df[!,"ActivePowerBalance__ACBus_sum"])-area2_df[!,"StateEstimationInjections__ACBus_sum"]
area3_df[!,"secheck"]=lag(area3_df[!,"ActivePowerBalance__ACBus_sum"])-area3_df[!,"StateEstimationInjections__ACBus_sum"]
println("area1 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("area2 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("area3 check max diff lag(ActivePowerBalance__ACBus_sum)-StateEstimationInjections__ACBus_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow_df=combine(groupby(bus_df, [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum)


[JuliaFormatter] reported by reviewdog 🐶

areaflow1_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="1"), areaflow_df)
areaflow2_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="2"), areaflow_df)
areaflow3_df=filter([:t,:area] => (t,area)-> (t>=1) && (area=="3"), areaflow_df)
areaflow1_df[!,"secheck"]=lag(areaflow1_df[!,"flowcontribution_sum"])-areaflow1_df[!,"loopflowcontribution_sum"]
areaflow2_df[!,"secheck"]=lag(areaflow2_df[!,"flowcontribution_sum"])-areaflow2_df[!,"loopflowcontribution_sum"]
areaflow3_df[!,"secheck"]=lag(areaflow3_df[!,"flowcontribution_sum"])-areaflow3_df[!,"loopflowcontribution_sum"]
println("areaflow1 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area1_df[!,"secheck"][2,:]))
println("areaflow2 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area2_df[!,"secheck"][2,:]))
println("areaflow3 check max diff lag(flowcontribution_sum)-loopflowcontribution_sum,", maximum(area3_df[!,"secheck"][2,:]))


[JuliaFormatter] reported by reviewdog 🐶

areaflow13_df=innerjoin(areaflow1_df, areaflow3_df, on=:t,renamecols = "_1" => "_3")
areaflow132_df=innerjoin(areaflow13_df, areaflow2_df, on=:t,renamecols = "" => "_2")
areaflow132_df[!,"flow_subsystem_a"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"loopflowcontribution_sum_2"]
areaflow132_df[!,"flow_subsystem_b"]=areaflow132_df[!,"loopflowcontribution_sum_1"]+areaflow132_df[!,"loopflowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
areaflow132_df[!,"seflow"]=areaflow132_df[!,"flowcontribution_sum_1"]+areaflow132_df[!,"flowcontribution_sum_3"]+areaflow132_df[!,"flowcontribution_sum_2"]
StatsPlots.@df(areaflow132_df,Plots.plot(:t,[:flow_subsystem_a :flow_subsystem_b :seflow]))


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,line,se=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()


[JuliaFormatter] reported by reviewdog 🐶

bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[], PowerLoad=[], ActivePowerBalance__ACBus=[],
StateEstimationInjections__ACBus=[], ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn], bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], ptdf[line,bn],
bus_ActivePowerBalance__ACBus[t,bn]*ptdf[line,bn],bus_StateEstimationInjections__ACBus[t,bn]*ptdf[line,bn]))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function check_models(uc0,uc2a)
sol0=Dict()
v0=uc0.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol0=Dict(zip(name.(v0),value.(v0)))


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end


[JuliaFormatter] reported by reviewdog 🐶

v01=uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")]


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end
v2a=uc2a.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol2a=Dict(zip(name.(v2a),value.(v2a)))


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in sol2a
if abs(sol2a[k]-sol0[k])>0.00001
println("k,",k,",",sol2a[k],",",sol0[k])
var0=variable_by_name(uc0.JuMPmodel,k)
fix(var0,sol2a[k],force=true)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶

end
return(objd0,sol0, objv0)
end


[JuliaFormatter] reported by reviewdog 🐶

function write_param_acbus(name,uc2a,line)
parameter_ACBus_2a=uc2a.parameters[InfrastructureSystems.Optimization.ParameterKey{PowerSimulationsDecomposition.StateEstimationInjections, ACBus}("")]
i=axes(parameter_ACBus_2a.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_ACBus_2a.parameter_array[j1,j2]
println(name,",j1,j2,",j1,",",j2,",value,",value(parameter_ACBus_2a.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_coeff(uc2a,linename,t,fnm,use_monitoredline)
fixedmw=0; fixedflow=0; areamw=Dict(); areaflow=Dict(); areaload=Dict(); arealoadflow=Dict()
for j=1:3 areamw[j]=0; areaflow[j]=0 end
c=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,t]


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(c[i], k)!=0


[JuliaFormatter] reported by reviewdog 🐶

if name(k)=="param"
areaload[i]=-normalized_coefficient(c[i], k)
println("area,",i,",",areaload[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if use_monitoredline==1
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][linename,t]
else
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, Line}("")][linename,t]
end
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
if name(k)=="param"
fixedmw=fixedmw+value(k);
fixedflow=fixedflow+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶

str=split(name(k),"{");
println("name,",name(k),",str,",str)
if first(str[1],1)=="A"
a= parse(Int64,first(str[2],1));# tt=firstindex(split(str[2],", "))
areamw[a]=areamw[a]+value(k)
areaflow[a]=areaflow[a]+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

return(fixedmw,areamw,fixedflow,areaflow)


[JuliaFormatter] reported by reviewdog 🐶

function write_param_internalload(name,uc0,line)
parameter_load_0=uc0.parameters[InfrastructureSystems.Optimization.ParameterKey{ActivePowerTimeSeriesParameter, PowerLoad}("")]
i=axes(parameter_load_0.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc0.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_load_0.parameter_array[j1,j2]
println(name,",load,j1,j2,",j1,",",j2,",value,",value(parameter_load_0.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end
function write_acbus(ucname,uc0,ax,t, debug=0)
area_param=Dict(); area_gen=Dict()
area_param[1]=0;area_param[2]=0;area_param[3]=0
area_gen[1]=0;area_gen[2]=0;area_gen[3]=0
for bus in ax[1]
for k in keys(PSI.get_expression(uc0, ActivePowerBalance(), PSY.ACBus)[bus,t].terms)
coef=PSI.get_expression(uc0, ActivePowerBalance(),PSY.ACBus)[bus,t].terms[k]
v=value(k);
if debug==1
println(ucname,";bus;t;",bus,";",t,";",k,";coefficient;", coef,";solution;",v,";coef*v;",coef*v)
end
n=Int64(floor(bus/100))
if name(k)=="param" area_param[n]=area_param[n]+coef*v else area_gen[n]=area_gen[n]+coef*v end


[JuliaFormatter] reported by reviewdog 🐶

end
return(area_gen,area_param)


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(get_bus(b))) == enforced_region[1] || get_name(get_area(get_bus(b))) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(b)) == enforced_region[1] || get_name(get_area(b)) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, "b", l)
end
end

Copy link

@github-actions github-actions bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function check_models(uc0,uc2a)
sol0=Dict()
v0=uc0.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol0=Dict(zip(name.(v0),value.(v0)))


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end


[JuliaFormatter] reported by reviewdog 🐶

v01=uc0.variables[InfrastructureSystems.Optimization.VariableKey{FlowActivePowerVariable, AreaInterchange}("")]


[JuliaFormatter] reported by reviewdog 🐶

fix(v,value.(v),force=true)
end
v2a=uc2a.variables[PowerSimulations.VariableKey{ActivePowerVariable, ThermalStandard}("")]
sol2a=Dict(zip(name.(v2a),value.(v2a)))


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in sol2a
if abs(sol2a[k]-sol0[k])>0.00001
println("k,",k,",",sol2a[k],",",sol0[k])
var0=variable_by_name(uc0.JuMPmodel,k)
fix(var0,sol2a[k],force=true)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶

end
return(objd0,sol0, objv0)
end


[JuliaFormatter] reported by reviewdog 🐶

function write_param_acbus(name,uc2a,line)
parameter_ACBus_2a=uc2a.parameters[InfrastructureSystems.Optimization.ParameterKey{PowerSimulationsDecomposition.StateEstimationInjections, ACBus}("")]
i=axes(parameter_ACBus_2a.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_ACBus_2a.parameter_array[j1,j2]
println(name,",j1,j2,",j1,",",j2,",value,",value(parameter_ACBus_2a.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_coeff(uc2a,linename,t,fnm,use_monitoredline)
fixedmw=0; fixedflow=0; areamw=Dict(); areaflow=Dict(); areaload=Dict(); arealoadflow=Dict()
for j=1:3 areamw[j]=0; areaflow[j]=0 end
c=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{CopperPlateBalanceConstraint, Area}("")][:,t]


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(c[i], k)!=0


[JuliaFormatter] reported by reviewdog 🐶

if name(k)=="param"
areaload[i]=-normalized_coefficient(c[i], k)
println("area,",i,",",areaload[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if use_monitoredline==1
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][linename,t]
else
con2a=uc2a.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, Line}("")][linename,t]
end
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
if name(k)=="param"
fixedmw=fixedmw+value(k);
fixedflow=fixedflow+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶

str=split(name(k),"{");
println("name,",name(k),",str,",str)
if first(str[1],1)=="A"
a= parse(Int64,first(str[2],1));# tt=firstindex(split(str[2],", "))
areamw[a]=areamw[a]+value(k)
areaflow[a]=areaflow[a]+value(k)*normalized_coefficient(con2a, k)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

return(fixedmw,areamw,fixedflow,areaflow)


[JuliaFormatter] reported by reviewdog 🐶

function write_param_internalload(name,uc0,line)
parameter_load_0=uc0.parameters[InfrastructureSystems.Optimization.ParameterKey{ActivePowerTimeSeriesParameter, PowerLoad}("")]
i=axes(parameter_load_0.parameter_array)


[JuliaFormatter] reported by reviewdog 🐶

for j1 in i[1]
con=uc0.constraints[InfrastructureSystems.Optimization.ConstraintKey{NetworkFlowConstraint, MonitoredLine}("")][line,j2]
param=parameter_load_0.parameter_array[j1,j2]
println(name,",load,j1,j2,",j1,",",j2,",value,",value(parameter_load_0.parameter_array[j1,j2]),",coeff,",normalized_coefficient(con,param))
end
end
end
function write_acbus(ucname,uc0,ax,t, debug=0)
area_param=Dict(); area_gen=Dict()
area_param[1]=0;area_param[2]=0;area_param[3]=0
area_gen[1]=0;area_gen[2]=0;area_gen[3]=0
for bus in ax[1]
for k in keys(PSI.get_expression(uc0, ActivePowerBalance(), PSY.ACBus)[bus,t].terms)
coef=PSI.get_expression(uc0, ActivePowerBalance(),PSY.ACBus)[bus,t].terms[k]
v=value(k);
if debug==1
println(ucname,";bus;t;",bus,";",t,";",k,";coefficient;", coef,";solution;",v,";coef*v;",coef*v)
end
n=Int64(floor(bus/100))
if name(k)=="param" area_param[n]=area_param[n]+coef*v else area_gen[n]=area_gen[n]+coef*v end


[JuliaFormatter] reported by reviewdog 🐶

end
return(area_gen,area_param)


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(get_bus(b))) == enforced_region[1] || get_name(get_area(get_bus(b))) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

if get_name(get_area(b)) == enforced_region[1] || get_name(get_area(b)) == enforced_region[2]


[JuliaFormatter] reported by reviewdog 🐶

for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, "b", l)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_cons_coeff(uc2a,con2a,fnm)
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_model(uc0,fnm)
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_constraints(uc0.JuMPmodel,; include_variable_in_set_constraints = true)
println(name(k),",",k)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,line="",se=0,standardload=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end
if standardload==1
ald = read_realized_variable(results_rt, "ActivePowerTimeSeriesParameter__StandardLoad")


[JuliaFormatter] reported by reviewdog 🐶

ald = read_realized_variable(results_rt, "ActivePowerTimeSeriesParameter__PowerLoad")
end
if line!="" am = read_realized_variable(results_rt, "FlowActivePowerVariable__MonitoredLine") end
ahvdc=[]
try
ahvdc=read_realized_variable(results_rt, "FlowActivePowerVariable__TwoTerminalHVDCLine")
catch e end
bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerVariable__RenewableDispatch=Dict()
bus_ActivePowerVariable__HydroDispatch=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()
bus_ActivePowerVariable__HVDC=Dict()
bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[],Renewable=[],Hydro=[],
PowerLoad=[], HVDC=[],ActivePowerBalance__ACBus=[], StateEstimationInjections__ACBus=[],
ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerVariable__RenewableDispatch[t,bn]=0
bus_ActivePowerVariable__HydroDispatch[t,bn]=0
bus_ActivePowerVariable__HVDC[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if standardload==1
try b=get_number(get_bus(get_component(StandardLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end
end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_from(get_arc(get_component(TwoTerminalHVDCLine, sys, i))))
for t=1:NT
#bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]-ahvdc[t,i]
bus_ActivePowerVariable__HVDC[t,b]=bus_ActivePowerVariable__HVDC[t,b]-ahvdc[t,i]
end
catch e println("error ",e,",i=",i) end
try b=get_number(get_to(get_arc(get_component(TwoTerminalHVDCLine, sys, i))))
for t=1:NT


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__HVDC[t,b]=bus_ActivePowerVariable__HVDC[t,b]+ahvdc[t,i]
end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(RenewableDispatch, sys, i)))
for t=1:NT bus_ActivePowerVariable__RenewableDispatch[t,b]=bus_ActivePowerVariable__RenewableDispatch[t,b]+ac1[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(HydroDispatch, sys, i)))
for t=1:NT bus_ActivePowerVariable__HydroDispatch[t,b]=bus_ActivePowerVariable__HydroDispatch[t,b]+ac2[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
if line=="" gsf=0 else gsf=ptdf[line,bn] end
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn],
bus_ActivePowerVariable__RenewableDispatch[t,bn],bus_ActivePowerVariable__HydroDispatch[t,bn],
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],bus_ActivePowerVariable__HVDC[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], gsf,
bus_ActivePowerBalance__ACBus[t,bn]*gsf,bus_StateEstimationInjections__ACBus[t,bn]*gsf))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"Renewable"]/100-bus_df[!,"Hydro"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

region=area_region_dict[get_name(get_area(get_bus(b)))]


[JuliaFormatter] reported by reviewdog 🐶

region=area_region_dict[get_name(get_area(b))]


[JuliaFormatter] reported by reviewdog 🐶

if length(selected_line)>0
for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, v, l)
end
end
end
if length(monitoredline_subsys) >0
for (k,s) in monitoredline_subsys
l=PSY.get_component(ACBranch, sys,k)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, v, l)
end
end


[JuliaFormatter] reported by reviewdog 🐶

tbus=get_to(get_arc(dc)); fbus=get_from(get_arc(dc))
tregion=area_region_dict[get_name(get_area(tbus))]
fregion=area_region_dict[get_name(get_area(fbus))]


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, tregion, dc)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

selected_line=["CA-1", "CB-1", "AB1"]
ln="A28"
limit=1.0 #2.0


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "2"),


[JuliaFormatter] reported by reviewdog 🐶

name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "2"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

area_dict=Dict("1" => "1", "2" => "2", "3" =>"3")
CF="NS3-0"
NS=3
netmodels=["AreaPTDFPowerModel","AreaPTDFPowerModel","AreaPTDFPowerModel"]
copperplates=[1,0,0]
NTS=[24,24,1]
area_region_dicts=[Dict(),
Dict(),
Dict()]
CF="NS3-1"
NS=3
netmodels=["AreaPTDFPowerModel","SplitAreaPTDFPowerModel","SplitAreaPTDFPowerModel"]
copperplates=[1,0,0]
NTS=[24,24,1]
area_region_dicts=[Dict(),
Dict("1" => "a", "2" => "b", "3" =>"c"),
Dict("1" => "a", "2" => "b", "3" =>"c")]
monitoredline_subsys=[Dict(),Dict(ln => ["a"]),Dict(ln => ["a"])]


[JuliaFormatter] reported by reviewdog 🐶

syss=[sys]
templates=[]
ptdf = VirtualPTDF(sys; tol = 1e-4, max_cache_size = 10000)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

NT=NTS[i]
if i>1
sys2=deepcopy(sys)
push!(syss,sys2)


[JuliaFormatter] reported by reviewdog 🐶

area_region_dict=area_region_dicts[i]
if netmodels[i]=="SplitAreaPTDFPowerModel"
subsys=unique(values(area_region_dict))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

template =
MultiProblemTemplate(NetworkModel(SplitAreaPTDFPowerModel; use_slacks=true,PTDF_matrix = ptdf,), subsys)
else
template = ProblemTemplate(NetworkModel(AreaPTDFPowerModel; use_slacks = true, PTDF_matrix = ptdf,))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if standardload==1
PSI.set_device_model!(template, StandardLoad , StaticPowerLoad)


[JuliaFormatter] reported by reviewdog 🐶

PSI.set_device_model!(template, DeviceModel(HydroDispatch, HPS.HydroDispatchReservoirBudget,
time_series_names = Dict{Any, String}(
PSI.ActivePowerTimeSeriesParameter => "max_active_power",
HPS.EnergyBudgetTimeSeriesParameter => "hydro_budget",)
)
)


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template, DeviceModel(MonitoredLine, StaticBranchUnbounded, use_slacks = true))
if NT>1


[JuliaFormatter] reported by reviewdog 🐶

end
if copperplates[i]==0
set_device_model!(template, DeviceModel(MonitoredLine, StaticBranch, use_slacks = true))


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if length(subsys)>1
#buildsubsystem1(sys2, area_region_dict, union(selected_line,monitoredlined_line))
buildsubsystem1(sys2, area_region_dict, selected_line,monitoredline_subsys[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if netmodels[i]=="SplitAreaPTDFPowerModel"
push!(dm, DecisionModel(MultiRegionProblem, templates[i], syss[i], name="UC"*string(i),optimizer=optimizer,
store_variable_names=true,initialize_model=false,optimizer_solve_log_print=false,
direct_mode_optimizer=true,check_numerical_bounds=false,
calculate_conflict=true,rebuild_model=false,system_to_file = true))
else
push!(dm, DecisionModel(templates[i],syss[i], name="UC"*string(i), optimizer=optimizer,
store_variable_names=true,initialize_model=false,optimizer_solve_log_print=false,
direct_mode_optimizer=true,check_numerical_bounds=false,
calculate_conflict=true,rebuild_model=false, system_to_file = true))
end
end
if NS==2
models = PSI.SimulationModels(decision_models=[dm[1],dm[2]])
else
models = PSI.SimulationModels(decision_models=[dm[1],dm[2],dm[3]])
end
uc_simulation_ffs=[]
feedforwards_dict=Dict()
for i in 1:NS-1


[JuliaFormatter] reported by reviewdog 🐶

#feed forward area interchange
FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

component_type=ThermalStandard, source=OnVariable, affected_values=[ActivePowerVariable],)
push!(uc_simulation_ff, SCFF)


[JuliaFormatter] reported by reviewdog 🐶

push!(uc_simulation_ffs,uc_simulation_ff)
feedforwards_dict["UC"*string(i+1)]=uc_simulation_ff


[JuliaFormatter] reported by reviewdog 🐶

results = SimulationResults(sim,ignore_status = true)


[JuliaFormatter] reported by reviewdog 🐶

results_uc=[]; uc_areainterchange=[];uc_monitoredline=[];
bus_df_uc=[];area_df=[];areaflow_df=[]
tm=Array(1:NTS[1])
af=DataFrames.DataFrame(;tm)


[JuliaFormatter] reported by reviewdog 🐶

push!(results_uc,get_decision_problem_results(results, "UC"*string(i)))
push!(uc_areainterchange,read_realized_variable(results_uc[i], "FlowActivePowerVariable__AreaInterchange"))
push!(uc_monitoredline,read_realized_variable(results_uc[i], "FlowActivePowerVariable__MonitoredLine"))
if netmodels[i]=="SplitAreaPTDFPowerModel"
push!(bus_df_uc,read_bus_df(results_uc[i], ln,1,standardload))
else
push!(bus_df_uc,read_bus_df(results_uc[i], ln,0,standardload))
end
bus_df_uc[i][!,"buscheck"]=bus_df_uc[i][!,"ActivePowerBalance__ACBus"]-bus_df_uc[i][!,"HVDC"]/100-bus_df_uc[i][!,"ThermalStandard"]/100-bus_df_uc[i][!,"Renewable"]/100-bus_df_uc[i][!,"Hydro"]/100-bus_df_uc[i][!,"PowerLoad"]/100
println("UC0 bus check ActivePowerBalance__ACBus<>Gen/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df_uc[i]))


[JuliaFormatter] reported by reviewdog 🐶

###### area MW check #########
#area level net Gen, load and ActivePowerBalance__ACBus_sum check
push!(area_df,combine(groupby(bus_df_uc[i], [:t, :area]), [:ThermalStandard, :Renewable, :Hydro, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum))
area_df[i][!,"areacheck"]=area_df[i][!,"ActivePowerBalance__ACBus_sum"]-area_df[i][!,"ThermalStandard_sum"]/100-area_df[i][!,"Renewable_sum"]/100-area_df[i][!,"Hydro_sum"]/100-area_df[i][!,"PowerLoad_sum"]/100
println("i,area_df check,",sum(area_df[i][!,"areacheck"])) #should be close to 0
area_df[i][!,"TotalGen"]=area_df[i][!,"ThermalStandard_sum"]+area_df[i][!,"Renewable_sum"]+area_df[i][!,"Hydro_sum"]
area_df[i][!,"NetInterchange"]=area_df[i][!,"TotalGen"]+area_df[i][!,"PowerLoad_sum"]
area_df[i][!,"InterchangeCheck"]=area_df[i][!,"NetInterchange"]-100*area_df[i][!,"ActivePowerBalance__ACBus_sum"]
println(area_df[i][:,["t","area","ActivePowerBalance__ACBus_sum","NetInterchange","InterchangeCheck"]])


[JuliaFormatter] reported by reviewdog 🐶

####### area flow check #########
#area flow contribution by area and interval
push!(areaflow_df, combine(groupby(bus_df_uc[i], [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

area_region_dict=area_dict
else
area_region_dict=area_region_dicts[i]
end
subsys=unique(values(area_region_dict))
calculated_flow_subsys=Dict()
calculated_loopflow_subsys=Dict()
NT=NTS[1]


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow_subsys_"*s]=zeros(NT)
af[!,"uc"*string(i)*"_loopflow_subsys_"*s]=zeros(NT)
end


[JuliaFormatter] reported by reviewdog 🐶

if netmodels[i]!="SplitAreaPTDFPowerModel"


[JuliaFormatter] reported by reviewdog 🐶

s=area_region_dict[r]
name="uc"*string(i)*"_flow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:flowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow"]= zeros(NT) #calculated_flow_uc0


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow"]=af[!,"uc"*string(i)*"_flow"]+af[!,"uc"*string(i)*"_flow_subsys_"*s]#calculated_flow_uc0_subsys[k]
end


[JuliaFormatter] reported by reviewdog 🐶

s=area_region_dict[r]
name="uc"*string(i)*"_flow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:flowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_loopflow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:loopflowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_seflow"]=zeros(NT)


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow_"*s]=zeros(NT)
af[!,"uc"*string(i)*"_seflow"]=af[!,"uc"*string(i)*"_seflow"]+af[!,"uc"*string(i)*"_flow_subsys_"*s]#calculated_flow_subsys[s]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_flow_"*s
af[!,name]=af[!,name]+af[!,"uc"*string(i)*"_flow_subsys_"*s1]#calculated_flow_subsys[s1]


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_flow_"*s
af[!,name]=af[!,name]+af[!,"uc"*string(i)*"_loopflow_subsys_"*s1]#calculated_loopflow_subsys[s1]
end
end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This files shouldn't be in the repositories

Copy link

@github-actions github-actions bot left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Remaining comments which cannot be posted as a review comment to avoid GitHub Rate Limit

JuliaFormatter

[JuliaFormatter] reported by reviewdog 🐶

function write_cons_coeff(uc2a,con2a,fnm)
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_variables(uc2a.JuMPmodel)
if normalized_coefficient(con2a, k)!=0
println("coef;",name(k),";",normalized_coefficient(con2a, k), ";is_fixed;",is_fixed(k),";value;",value(k))
end
end


[JuliaFormatter] reported by reviewdog 🐶

function write_model(uc0,fnm)
open(fnm,"w") do io


[JuliaFormatter] reported by reviewdog 🐶

for k in all_constraints(uc0.JuMPmodel,; include_variable_in_set_constraints = true)
println(name(k),",",k)
end
end


[JuliaFormatter] reported by reviewdog 🐶

function objterms(uc0,ff)
objd0=Dict(); sol0=Dict(); objv0=0
obj0=objective_function(uc0.JuMPmodel); objd0=Dict(zip(name.(keys(obj0.terms)),values(obj0.terms)));
av0 = JuMP.all_variables(uc0.JuMPmodel); sol0=Dict(zip(name.(av0), value.(av0)))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for (k,v) in objd0
objv0=objv0+sol0[k]*objd0[k]
println("name;value;objcoef;obj_contribution;",k,";",sol0[k],";",objd0[k],";", sol0[k]*objd0[k])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

function read_bus_df(results_rt,line="",se=0,standardload=0)


[JuliaFormatter] reported by reviewdog 🐶

if se==1 al = read_realized_variable(results_rt, "StateEstimationInjections__ACBus") end
if standardload==1
ald = read_realized_variable(results_rt, "ActivePowerTimeSeriesParameter__StandardLoad")


[JuliaFormatter] reported by reviewdog 🐶

ald = read_realized_variable(results_rt, "ActivePowerTimeSeriesParameter__PowerLoad")
end
if line!="" am = read_realized_variable(results_rt, "FlowActivePowerVariable__MonitoredLine") end
ahvdc=[]
try
ahvdc=read_realized_variable(results_rt, "FlowActivePowerVariable__TwoTerminalHVDCLine")
catch e end
bus_ActivePowerVariable__ThermalStandard=Dict()
bus_ActivePowerVariable__RenewableDispatch=Dict()
bus_ActivePowerVariable__HydroDispatch=Dict()
bus_ActivePowerTimeSeriesParameter__PowerLoad=Dict()
bus_ActivePowerBalance__ACBus=Dict()
bus_StateEstimationInjections__ACBus=Dict()
bus_ActivePowerVariable__HVDC=Dict()
bus_df=DataFrames.DataFrame(;
t =Int64[], bus=[], area = [], ThermalStandard=[],Renewable=[],Hydro=[],
PowerLoad=[], HVDC=[],ActivePowerBalance__ACBus=[], StateEstimationInjections__ACBus=[],
ptdf=[], flowcontribution=[], loopflowcontribution=[]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); bus_ActivePowerVariable__ThermalStandard[t,bn]=0
bus_ActivePowerVariable__RenewableDispatch[t,bn]=0
bus_ActivePowerVariable__HydroDispatch[t,bn]=0
bus_ActivePowerVariable__HVDC[t,bn]=0
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn]=0; bus_ActivePowerBalance__ACBus[t,bn]=0
bus_StateEstimationInjections__ACBus[t,bn]=0


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_ActivePowerBalance__ACBus[t,b]=ab[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=parse(Int,i)
for t=1:NT bus_StateEstimationInjections__ACBus[t,b]=al[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if standardload==1
try b=get_number(get_bus(get_component(StandardLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(PowerLoad, sys, i)))
for t=1:NT bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]=bus_ActivePowerTimeSeriesParameter__PowerLoad[t,b]+ald[t,i] end
catch e println("error= ",e,",i=",i) end
end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_from(get_arc(get_component(TwoTerminalHVDCLine, sys, i))))
for t=1:NT
#bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]-ahvdc[t,i]
bus_ActivePowerVariable__HVDC[t,b]=bus_ActivePowerVariable__HVDC[t,b]-ahvdc[t,i]
end
catch e println("error ",e,",i=",i) end
try b=get_number(get_to(get_arc(get_component(TwoTerminalHVDCLine, sys, i))))
for t=1:NT


[JuliaFormatter] reported by reviewdog 🐶

bus_ActivePowerVariable__HVDC[t,b]=bus_ActivePowerVariable__HVDC[t,b]+ahvdc[t,i]
end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(ThermalStandard, sys, i)))
for t=1:NT bus_ActivePowerVariable__ThermalStandard[t,b]=bus_ActivePowerVariable__ThermalStandard[t,b]+ac[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(RenewableDispatch, sys, i)))
for t=1:NT bus_ActivePowerVariable__RenewableDispatch[t,b]=bus_ActivePowerVariable__RenewableDispatch[t,b]+ac1[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

try b=get_number(get_bus(get_component(HydroDispatch, sys, i)))
for t=1:NT bus_ActivePowerVariable__HydroDispatch[t,b]=bus_ActivePowerVariable__HydroDispatch[t,b]+ac2[t,i] end
catch e println("error ",e,",i=",i) end


[JuliaFormatter] reported by reviewdog 🐶

for b in PSY.get_components(PSY.Bus,sys)
bn=get_number(b); area=get_name(get_area(b));
if line=="" gsf=0 else gsf=ptdf[line,bn] end
push!(bus_df,(t,bn,area, bus_ActivePowerVariable__ThermalStandard[t,bn],
bus_ActivePowerVariable__RenewableDispatch[t,bn],bus_ActivePowerVariable__HydroDispatch[t,bn],
bus_ActivePowerTimeSeriesParameter__PowerLoad[t,bn],bus_ActivePowerVariable__HVDC[t,bn],
bus_ActivePowerBalance__ACBus[t,bn], bus_StateEstimationInjections__ACBus[t,bn], gsf,
bus_ActivePowerBalance__ACBus[t,bn]*gsf,bus_StateEstimationInjections__ACBus[t,bn]*gsf))


[JuliaFormatter] reported by reviewdog 🐶

bus_df[!,"buscheck"]=bus_df[!,"ActivePowerBalance__ACBus"]-bus_df[!,"ThermalStandard"]/100-bus_df[!,"Renewable"]/100-bus_df[!,"Hydro"]/100-bus_df[!,"PowerLoad"]/100
println("bus check ActivePowerBalance__ACBus<>ThermalStandard/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

region=area_region_dict[get_name(get_area(get_bus(b)))]


[JuliaFormatter] reported by reviewdog 🐶

region=area_region_dict[get_name(get_area(b))]


[JuliaFormatter] reported by reviewdog 🐶

if length(selected_line)>0
for b in selected_line
l=PSY.get_component(ACBranch, sys,b)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, v, l)
end
end
end
if length(monitoredline_subsys) >0
for (k,s) in monitoredline_subsys
l=PSY.get_component(ACBranch, sys,k)


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, v, l)
end
end


[JuliaFormatter] reported by reviewdog 🐶

tbus=get_to(get_arc(dc)); fbus=get_from(get_arc(dc))
tregion=area_region_dict[get_name(get_area(tbus))]
fregion=area_region_dict[get_name(get_area(fbus))]


[JuliaFormatter] reported by reviewdog 🐶

add_component_to_subsystem!(sys, tregion, dc)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

selected_line=["CA-1", "CB-1", "AB1"]
ln="A28"
limit=1.0 #2.0


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

name="1_2", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "2"),


[JuliaFormatter] reported by reviewdog 🐶

name="1_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "1"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

name="2_3", available=true, active_power_flow=0.0, from_area=get_component(Area, sys, "2"), to_area=get_component(Area, sys, "3"),


[JuliaFormatter] reported by reviewdog 🐶

area_dict=Dict("1" => "1", "2" => "2", "3" =>"3")
CF="NS3-0"
NS=3
netmodels=["AreaPTDFPowerModel","AreaPTDFPowerModel","AreaPTDFPowerModel"]
copperplates=[1,0,0]
NTS=[24,24,1]
area_region_dicts=[Dict(),
Dict(),
Dict()]
CF="NS3-1"
NS=3
netmodels=["AreaPTDFPowerModel","SplitAreaPTDFPowerModel","SplitAreaPTDFPowerModel"]
copperplates=[1,0,0]
NTS=[24,24,1]
area_region_dicts=[Dict(),
Dict("1" => "a", "2" => "b", "3" =>"c"),
Dict("1" => "a", "2" => "b", "3" =>"c")]
monitoredline_subsys=[Dict(),Dict(ln => ["a"]),Dict(ln => ["a"])]


[JuliaFormatter] reported by reviewdog 🐶

syss=[sys]
templates=[]
ptdf = VirtualPTDF(sys; tol = 1e-4, max_cache_size = 10000)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

NT=NTS[i]
if i>1
sys2=deepcopy(sys)
push!(syss,sys2)


[JuliaFormatter] reported by reviewdog 🐶

area_region_dict=area_region_dicts[i]
if netmodels[i]=="SplitAreaPTDFPowerModel"
subsys=unique(values(area_region_dict))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

template =
MultiProblemTemplate(NetworkModel(SplitAreaPTDFPowerModel; use_slacks=true,PTDF_matrix = ptdf,), subsys)
else
template = ProblemTemplate(NetworkModel(AreaPTDFPowerModel; use_slacks = true, PTDF_matrix = ptdf,))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if standardload==1
PSI.set_device_model!(template, StandardLoad , StaticPowerLoad)


[JuliaFormatter] reported by reviewdog 🐶

PSI.set_device_model!(template, DeviceModel(HydroDispatch, HPS.HydroDispatchReservoirBudget,
time_series_names = Dict{Any, String}(
PSI.ActivePowerTimeSeriesParameter => "max_active_power",
HPS.EnergyBudgetTimeSeriesParameter => "hydro_budget",)
)
)


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template, DeviceModel(MonitoredLine, StaticBranchUnbounded, use_slacks = true))
if NT>1


[JuliaFormatter] reported by reviewdog 🐶

end
if copperplates[i]==0
set_device_model!(template, DeviceModel(MonitoredLine, StaticBranch, use_slacks = true))


[JuliaFormatter] reported by reviewdog 🐶

set_device_model!(template,
DeviceModel(Line, StaticBranchUnbounded; attributes=Dict("filter_function" => x -> get_name(x) in selected_line),))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if length(subsys)>1
#buildsubsystem1(sys2, area_region_dict, union(selected_line,monitoredlined_line))
buildsubsystem1(sys2, area_region_dict, selected_line,monitoredline_subsys[i])


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

if netmodels[i]=="SplitAreaPTDFPowerModel"
push!(dm, DecisionModel(MultiRegionProblem, templates[i], syss[i], name="UC"*string(i),optimizer=optimizer,
store_variable_names=true,initialize_model=false,optimizer_solve_log_print=false,
direct_mode_optimizer=true,check_numerical_bounds=false,
calculate_conflict=true,rebuild_model=false,system_to_file = true))
else
push!(dm, DecisionModel(templates[i],syss[i], name="UC"*string(i), optimizer=optimizer,
store_variable_names=true,initialize_model=false,optimizer_solve_log_print=false,
direct_mode_optimizer=true,check_numerical_bounds=false,
calculate_conflict=true,rebuild_model=false, system_to_file = true))
end
end
if NS==2
models = PSI.SimulationModels(decision_models=[dm[1],dm[2]])
else
models = PSI.SimulationModels(decision_models=[dm[1],dm[2],dm[3]])
end
uc_simulation_ffs=[]
feedforwards_dict=Dict()
for i in 1:NS-1


[JuliaFormatter] reported by reviewdog 🐶

#feed forward area interchange
FVFF_area_interchange = FixValueFeedforward(;component_type=AreaInterchange,source=FlowActivePowerVariable,affected_values=[FlowActivePowerVariable],)


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

component_type=ThermalStandard, source=OnVariable, affected_values=[ActivePowerVariable],)
push!(uc_simulation_ff, SCFF)


[JuliaFormatter] reported by reviewdog 🐶

push!(uc_simulation_ffs,uc_simulation_ff)
feedforwards_dict["UC"*string(i+1)]=uc_simulation_ff


[JuliaFormatter] reported by reviewdog 🐶

results = SimulationResults(sim,ignore_status = true)


[JuliaFormatter] reported by reviewdog 🐶

results_uc=[]; uc_areainterchange=[];uc_monitoredline=[];
bus_df_uc=[];area_df=[];areaflow_df=[]
tm=Array(1:NTS[1])
af=DataFrames.DataFrame(;tm)


[JuliaFormatter] reported by reviewdog 🐶

push!(results_uc,get_decision_problem_results(results, "UC"*string(i)))
push!(uc_areainterchange,read_realized_variable(results_uc[i], "FlowActivePowerVariable__AreaInterchange"))
push!(uc_monitoredline,read_realized_variable(results_uc[i], "FlowActivePowerVariable__MonitoredLine"))
if netmodels[i]=="SplitAreaPTDFPowerModel"
push!(bus_df_uc,read_bus_df(results_uc[i], ln,1,standardload))
else
push!(bus_df_uc,read_bus_df(results_uc[i], ln,0,standardload))
end
bus_df_uc[i][!,"buscheck"]=bus_df_uc[i][!,"ActivePowerBalance__ACBus"]-bus_df_uc[i][!,"HVDC"]/100-bus_df_uc[i][!,"ThermalStandard"]/100-bus_df_uc[i][!,"Renewable"]/100-bus_df_uc[i][!,"Hydro"]/100-bus_df_uc[i][!,"PowerLoad"]/100
println("UC0 bus check ActivePowerBalance__ACBus<>Gen/100-PowerLoad/100 ,",filter([:t,:buscheck] => (t,buscheck)-> (t>=1) && (abs(buscheck)>0.000001), bus_df_uc[i]))


[JuliaFormatter] reported by reviewdog 🐶

###### area MW check #########
#area level net Gen, load and ActivePowerBalance__ACBus_sum check
push!(area_df,combine(groupby(bus_df_uc[i], [:t, :area]), [:ThermalStandard, :Renewable, :Hydro, :PowerLoad, :ActivePowerBalance__ACBus, :StateEstimationInjections__ACBus] .=> sum))
area_df[i][!,"areacheck"]=area_df[i][!,"ActivePowerBalance__ACBus_sum"]-area_df[i][!,"ThermalStandard_sum"]/100-area_df[i][!,"Renewable_sum"]/100-area_df[i][!,"Hydro_sum"]/100-area_df[i][!,"PowerLoad_sum"]/100
println("i,area_df check,",sum(area_df[i][!,"areacheck"])) #should be close to 0
area_df[i][!,"TotalGen"]=area_df[i][!,"ThermalStandard_sum"]+area_df[i][!,"Renewable_sum"]+area_df[i][!,"Hydro_sum"]
area_df[i][!,"NetInterchange"]=area_df[i][!,"TotalGen"]+area_df[i][!,"PowerLoad_sum"]
area_df[i][!,"InterchangeCheck"]=area_df[i][!,"NetInterchange"]-100*area_df[i][!,"ActivePowerBalance__ACBus_sum"]
println(area_df[i][:,["t","area","ActivePowerBalance__ACBus_sum","NetInterchange","InterchangeCheck"]])


[JuliaFormatter] reported by reviewdog 🐶

####### area flow check #########
#area flow contribution by area and interval
push!(areaflow_df, combine(groupby(bus_df_uc[i], [:t, :area]), [:flowcontribution, :loopflowcontribution] .=> sum))


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

area_region_dict=area_dict
else
area_region_dict=area_region_dicts[i]
end
subsys=unique(values(area_region_dict))
calculated_flow_subsys=Dict()
calculated_loopflow_subsys=Dict()
NT=NTS[1]


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow_subsys_"*s]=zeros(NT)
af[!,"uc"*string(i)*"_loopflow_subsys_"*s]=zeros(NT)
end


[JuliaFormatter] reported by reviewdog 🐶

if netmodels[i]!="SplitAreaPTDFPowerModel"


[JuliaFormatter] reported by reviewdog 🐶

s=area_region_dict[r]
name="uc"*string(i)*"_flow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:flowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow"]= zeros(NT) #calculated_flow_uc0


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow"]=af[!,"uc"*string(i)*"_flow"]+af[!,"uc"*string(i)*"_flow_subsys_"*s]#calculated_flow_uc0_subsys[k]
end


[JuliaFormatter] reported by reviewdog 🐶

s=area_region_dict[r]
name="uc"*string(i)*"_flow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:flowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_loopflow_subsys_"*s
af[!,name]=af[!,name]+filter([:t,:area]=>(t,area)->(area==r),areaflow_df[i])[!,:loopflowcontribution_sum]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_seflow"]=zeros(NT)


[JuliaFormatter] reported by reviewdog 🐶

af[!,"uc"*string(i)*"_flow_"*s]=zeros(NT)
af[!,"uc"*string(i)*"_seflow"]=af[!,"uc"*string(i)*"_seflow"]+af[!,"uc"*string(i)*"_flow_subsys_"*s]#calculated_flow_subsys[s]


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_flow_"*s
af[!,name]=af[!,name]+af[!,"uc"*string(i)*"_flow_subsys_"*s1]#calculated_flow_subsys[s1]


[JuliaFormatter] reported by reviewdog 🐶

name="uc"*string(i)*"_flow_"*s
af[!,name]=af[!,name]+af[!,"uc"*string(i)*"_loopflow_subsys_"*s1]#calculated_loopflow_subsys[s1]
end
end


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶


[JuliaFormatter] reported by reviewdog 🐶

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants