AC powerflow breaks in the presence of switches

[eperim]

It looks like AC powerflow does not work for networks that have switches, even for the examples contained within the project itself:

julia> network_data = PowerModels.parse_file("matpower/case5_sw.m")
Dict{String, Any} with 13 entries:
  "bus"            => Dict{String, Any}("4"=>Dict{String, Any}("zone"=>1, "bus_i"=>4, "bus_type"=>3, "vmax"=>1.1, "source_id"=>Any["bus", 4], "area"=>1, "vmin"=>0.9, "index"=>4, "va"=>0.0, "vm"=>1.0…), "1"=>Dict{String, Any}("zone"=>1, "bus_i"=>1, "bus_type"=>2, "vmax"=>1.1, "source_id"=>Any["bus", 1], "area"=>1, "vmin"=>0.9, "index"=>1, "va"=>0.048935, "vm"=>…
  "source_type"    => "matpower"
  "name"           => "case5"
  "dcline"         => Dict{String, Any}()
  "source_version" => "2"
  "gen"            => Dict{String, Any}("4"=>Dict{String, Any}("pg"=>0.0, "model"=>2, "shutdown"=>0.0, "startup"=>0.0, "qg"=>-0.1, "gen_bus"=>4, "pmax"=>2.0, "vg"=>1.0, "mbase"=>100.0, "source_id"=>Any["gen", 4]…), "1"=>Dict{String, Any}("pg"=>0.4, "model"=>2, "shutdown"=>0.0, "startup"=>0.0, "qg"=>0.3, "gen_bus"=>1, "pmax"=>0.4, "vg"=>1.0, "mbase"=>100.0, "so…
  "branch"         => Dict{String, Any}("4"=>Dict{String, Any}("br_r"=>0.00297, "rate_a"=>4.26, "shift"=>-0.0174533, "rate_b"=>4.26, "br_x"=>0.0297, "rate_c"=>4.26, "g_to"=>0.0, "g_fr"=>0.0, "source_id"=>Any["branch", 4], "b_fr"=>0.00337…), "1"=>Dict{String, Any}("br_r"=>0.00304, "rate_a"=>4.26, "shift"=>0.0, "rate_b"=>4.26, "br_x"=>0.0304, "rate_c"=>4.26, "g_…
  "storage"        => Dict{String, Any}()
  "switch"         => Dict{String, Any}("1"=>Dict{String, Any}("qsw"=>0.9861, "source_id"=>Any["switch", 1], "f_bus"=>1, "thermal_rating"=>10.0, "status"=>1, "t_bus"=>2, "psw"=>3.0, "index"=>1, "state"=>1), "2"=>Dict{String, Any}("qsw"=>0.0, "source_id"=>Any["switch", 2], "f_bus"=>3, "thermal_rating"=>10.0, "status"=>1, "t_bus"=>2, "psw"=>0.0, "index"=>2, "sta…
  "baseMVA"        => 100.0
  "per_unit"       => true
  "shunt"          => Dict{String, Any}("1"=>Dict{String, Any}("source_id"=>Any["bus", 2], "shunt_bus"=>2, "status"=>1, "gs"=>0.05, "bs"=>0.1, "index"=>1))
  "load"           => Dict{String, Any}("1"=>Dict{String, Any}("source_id"=>Any["bus", 2], "load_bus"=>2, "status"=>1, "qd"=>0.9861, "pd"=>3.0, "index"=>1), "2"=>Dict{String, Any}("source_id"=>Any["bus", 3], "load_bus"=>3, "status"=>1, "qd"=>0.9861, "pd"=>3.0, "index"=>2), "3"=>Dict{String, Any}("source_id"=>Any["bus", 4], "load_bus"=>4, "status"=>1, "qd"=>1.3…

julia> network_data["switch"]["1"]
Dict{String, Any} with 9 entries:
  "qsw"            => 0.9861
  "source_id"      => Any["switch", 1]
  "f_bus"          => 1
  "thermal_rating" => 10.0
  "status"         => 1
  "t_bus"          => 2
  "psw"            => 3.0
  "index"          => 1
  "state"          => 1

julia> sol = solve_ac_pf(network_data, Ipopt.Optimizer)
[error | PowerModels]: active power decision variables appear to be missing for switch components
ERROR: active power decision variables appear to be missing for switch components
Stacktrace:
  [1] error(logger::Memento.Logger, msg::String)
    @ Memento ~/.julia/packages/Memento/xnHxE/src/loggers.jl:442
  [2] _check_var_keys(vars::Dict{Any, Any}, keys::Vector{Tuple{Int64, Int64, Int64}}, var_name::String, comp_name::String)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/core/constraint.jl:10
  [3] constraint_power_balance(pm::ACPPowerModel, n::Int64, i::Int64, bus_arcs::Vector{Tuple{…}}, bus_arcs_dc::Vector{Tuple{…}}, bus_arcs_sw::Vector{Tuple{…}}, bus_gens::Vector{Int64}, bus_storage::Vector{Int64}, bus_pd::Dict{Int64, Float64}, bus_qd::Dict{Int64, Float64}, bus_gs::Dict{Int64, Float64}, bus_bs::Dict{Int64, Float64})
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/form/acp.jl:61
  [4] constraint_power_balance(pm::ACPPowerModel, i::Int64; nw::Int64)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/core/constraint_template.jl:190
  [5] constraint_power_balance
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/constraint_template.jl:174 [inlined]
  [6] build_pf(pm::ACPPowerModel)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:44
  [7] instantiate_model(data::Dict{String, Any}, model_type::Type, build_method::typeof(build_pf), ref_add_core!::typeof(ref_add_core!), global_keys::Set{String}, it::Symbol; ref_extensions::Vector{Any}, kwargs::@Kwargs{})
    @ InfrastructureModels ~/.julia/packages/InfrastructureModels/C2xBM/src/core/base.jl:370
  [8] instantiate_model
    @ ~/.julia/packages/InfrastructureModels/C2xBM/src/core/base.jl:350 [inlined]
  [9] #instantiate_model#299
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:52 [inlined]
 [10] instantiate_model
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:51 [inlined]
 [11] solve_model(data::Dict{String, Any}, model_type::Type, optimizer::Type, build_method::typeof(build_pf); ref_extensions::Vector{Any}, solution_processors::Vector{Any}, relax_integrality::Bool, multinetwork::Bool, kwargs::@Kwargs{})
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:33
 [12] solve_model
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:22 [inlined]
 [13] #solve_pf#984
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:13 [inlined]
 [14] solve_pf
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:12 [inlined]
 [15] #solve_ac_pf#982
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:3 [inlined]
 [16] solve_ac_pf(file::Dict{String, Any}, optimizer::Type)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:2
 [17] top-level scope
    @ REPL[23]:1

This also happens if I call set_ac_pf_start_values!(network_data) before calling solve_ac_pf. It looks to me that all relevant fields for the switch are there, so I'd imagine the issue is somewhere around the instantiation of the model? Or am I doing something wrong?

[ccoffrin]

You should call one of these OPF variants that support switch data, https://github.com/lanl-ansi/PowerModels.jl/blob/master/src/prob/test.jl#L49-L221

[eperim]

Is there any implementation of powerflow that supports switches?

[eperim]

@ccoffrin Just confirming: is there any method in the package that supports running PF with switches, or does PowerModels only support switches in OPF runs?

[ccoffrin]

Yes that is correct. I personally have not thought deeply about what the best power flow formulation would be in the presence of switch data. I expect they would introduce a bunch of numerical issues in the Jacobean that would need to be handled with some care. That said, I do think it would be quick and easy to use the OPF example as a guide to extend this PF formulation (https://github.com/lanl-ansi/PowerModels.jl/blob/master/src/prob/pf.jl#L16-L73) to support switch data.

One approach I can recommend is to call resolve_swithces! (spelled wrong! #936) on the data and then calling any of the standard power flow methods. This method uses the status of the switches in the dataset to convert the network into a standard bus-branch model.

https://github.com/lanl-ansi/PowerModels.jl/blob/master/src/core/data.jl#L2686-L2749

[frederikgeth]

fwiw, this approach is used for state estimation with switches, i.e. augment the voltage vector with p/q variables for each switch: https://ieeexplore.ieee.org/abstract/document/709101. In my understanding that avoids problem with the jacobian when modeling very small impedances.

[eperim]

Thanks a lot for the help. However, I am still having an issue when trying to use the proposed approach with resolve_swithces!:

julia> network_data = PowerModels.parse_file("matpower/case5_sw.m")
Dict{String, Any} with 13 entries:
  "bus"            => Dict{String, Any}("4"=>Dict{String, Any}("zone"=>1, "bus_i"=>4, "bus_type"=>3, "vmax"=>1.1, "source_id"=>Any["bus", 4], "area"=>1, "vmin"=>0.9, "index"=>4, "va"=>0.0, "vm"=>1.0…), "1"=>Dict{String, Any}("zone"=>1, "bus_i"=>1, "bus_type"=>2, "vmax"=>1.1, "source_id"=>Any["bus", 1], "area"=>1, "vmin"=>0.9, "index"=>1, "va"=>0.048935, "vm"=>…
  "source_type"    => "matpower"
  "name"           => "case5"
  "dcline"         => Dict{String, Any}()
  "source_version" => "2"
  "gen"            => Dict{String, Any}("4"=>Dict{String, Any}("pg"=>0.0, "model"=>2, "shutdown"=>0.0, "startup"=>0.0, "qg"=>-0.1, "gen_bus"=>4, "pmax"=>2.0, "vg"=>1.0, "mbase"=>100.0, "source_id"=>Any["gen", 4]…), "1"=>Dict{String, Any}("pg"=>0.4, "model"=>2, "shutdown"=>0.0, "startup"=>0.0, "qg"=>0.3, "gen_bus"=>1, "pmax"=>0.4, "vg"=>1.0, "mbase"=>100.0, "so…
  "branch"         => Dict{String, Any}("4"=>Dict{String, Any}("br_r"=>0.00297, "rate_a"=>4.26, "shift"=>-0.0174533, "rate_b"=>4.26, "br_x"=>0.0297, "rate_c"=>4.26, "g_to"=>0.0, "g_fr"=>0.0, "source_id"=>Any["branch", 4], "b_fr"=>0.00337…), "1"=>Dict{String, Any}("br_r"=>0.00304, "rate_a"=>4.26, "shift"=>0.0, "rate_b"=>4.26, "br_x"=>0.0304, "rate_c"=>4.26, "g_…
  "storage"        => Dict{String, Any}()
  "switch"         => Dict{String, Any}("1"=>Dict{String, Any}("qsw"=>0.9861, "source_id"=>Any["switch", 1], "f_bus"=>1, "thermal_rating"=>10.0, "status"=>1, "t_bus"=>2, "psw"=>3.0, "index"=>1, "state"=>1), "2"=>Dict{String, Any}("qsw"=>0.0, "source_id"=>Any["switch", 2], "f_bus"=>3, "thermal_rating"=>10.0, "status"=>1, "t_bus"=>2, "psw"=>0.0, "index"=>2, "sta…
  "baseMVA"        => 100.0
  "per_unit"       => true
  "shunt"          => Dict{String, Any}("1"=>Dict{String, Any}("source_id"=>Any["bus", 2], "shunt_bus"=>2, "status"=>1, "gs"=>0.05, "bs"=>0.1, "index"=>1))
  "load"           => Dict{String, Any}("1"=>Dict{String, Any}("source_id"=>Any["bus", 2], "load_bus"=>2, "status"=>1, "qd"=>0.9861, "pd"=>3.0, "index"=>1), "2"=>Dict{String, Any}("source_id"=>Any["bus", 3], "load_bus"=>3, "status"=>1, "qd"=>0.9861, "pd"=>3.0, "index"=>2), "3"=>Dict{String, Any}("source_id"=>Any["bus", 4], "load_bus"=>4, "status"=>1, "qd"=>1.3…

julia> resolve_swithces!(network_data)
[info | PowerModels]: merged buses 2,1 in to bus 1 based on switch status
[info | PowerModels]: removed 3 switch components
Dict{String, Any}()

julia> sol = solve_ac_pf(network_data, Ipopt.Optimizer)
ERROR: AssertionError: bus["bus_type"] == 2
Stacktrace:
  [1] build_pf(pm::ACPPowerModel)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:49
  [2] instantiate_model(data::Dict{String, Any}, model_type::Type, build_method::typeof(build_pf), ref_add_core!::typeof(ref_add_core!), global_keys::Set{String}, it::Symbol; ref_extensions::Vector{Any}, kwargs::@Kwargs{})
    @ InfrastructureModels ~/.julia/packages/InfrastructureModels/C2xBM/src/core/base.jl:370
  [3] instantiate_model
    @ ~/.julia/packages/InfrastructureModels/C2xBM/src/core/base.jl:350 [inlined]
  [4] #instantiate_model#299
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:52 [inlined]
  [5] instantiate_model
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:51 [inlined]
  [6] solve_model(data::Dict{String, Any}, model_type::Type, optimizer::Type, build_method::typeof(build_pf); ref_extensions::Vector{Any}, solution_processors::Vector{Any}, relax_integrality::Bool, multinetwork::Bool, kwargs::@Kwargs{})
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:33
  [7] solve_model
    @ ~/.julia/packages/PowerModels/Dhfhc/src/core/base.jl:22 [inlined]
  [8] #solve_pf#984
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:13 [inlined]
  [9] solve_pf
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:12 [inlined]
 [10] #solve_ac_pf#982
    @ ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:3 [inlined]
 [11] solve_ac_pf(file::Dict{String, Any}, optimizer::Type)
    @ PowerModels ~/.julia/packages/PowerModels/Dhfhc/src/prob/pf.jl:2
 [12] top-level scope
    @ REPL[50]:1

My impression is that when the buses get merged, the bus type is not updated accordingly. In this case, buses 1 and 2 get merged. 1 is PQ, while 2 is PV, so I reckon the result should be a PV bus, but since 1 is used as reference, the resulting bus is labelled as PQ, which breaks the assertion. I imagine that we would need to update the bus type accordingly in this loop. Do you agree with this?

[KSepetanc]

Just asking to check if the bus merge behavior is limited to PF computation/commands. I do OPF and want to preserve topology as much as possible. Does any of parse_file or build_ref do the merge?

[eperim]

Just asking to check if the bus merge behavior is limited to PF computation/commands. I do OPF and want to preserve topology as much as possible.

Take my opinion with a grain of salt, as I am just a user, but I believe this bus type issue only really matters for PF calculations, as PV and PQ buses are treated in intrinsically different manners during the PF iterations. That shouldn't be the case for OPF when it comes to PV/PQ buses (the reference bus will still be treated differently because its angle is set to zero).