Jd/ac power flow refactor

Project.toml

@@ -10,7 +10,6 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
 PowerNetworkMatrices = "bed98974-b02a-5e2f-9fe0-a103f5c450dd"
 PowerSystems = "bcd98974-b02a-5e2f-9ee0-a103f5c450dd"
-SharedArrays = "1a1011a3-84de-559e-8e89-a11a2f7dc383"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 
 [compat]

docs/src/modeler_guide/power_flow.md

@@ -67,7 +67,7 @@ For instance, initializing dynamic simulations. Also, because [`run_powerflow!`]
 are also available for [`run_powerflow!`](@ref)
 
 ````@example generated_power_flow
-run_powerflow!(system_data; finite_diff = true, method = :newton)
+run_powerflow!(system_data; method = :newton)
 ````
 
 After running the power flow command this are the values of the

src/PowerFlowData.jl

@@ -1,60 +1,61 @@
 """
-Structure containing all the data required for the evaluation of the power 
+Structure containing all the data required for the evaluation of the power
 flows and angles, as well as these ones.
 
 # Arguments:
 - `bus_lookup::Dict{Int, Int}`:
-        dictionary linking the system's bus number with the rows of either 
+        dictionary linking the system's bus number with the rows of either
         "power_network_matrix" or "aux_network_matrix".
 - `branch_lookup::Dict{String, Int}`:
-        dictionary linking the branch name with the column name of either the 
+        dictionary linking the branch name with the column name of either the
         "power_network_matrix" or "aux_network_matrix".
 - `bus_activepower_injection::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus active power injection. b: number of 
+        "(b, t)" matrix containing the bus active power injection. b: number of
         buses, t: number of time period.
 - `bus_reactivepower_injection::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus reactive power injection. b: number 
+        "(b, t)" matrix containing the bus reactive power injection. b: number
         of buses, t: number of time period.
 - `bus_activepower_withdrawals::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus reactive power withdrawals. b: 
+        "(b, t)" matrix containing the bus reactive power withdrawals. b:
         number of buses, t: number of time period.
 - `bus_reactivepower_withdrawals::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus reactive power withdrawals. b: 
+        "(b, t)" matrix containing the bus reactive power withdrawals. b:
         number of buses, t: number of time period.
+- `bus_reactivepower_bounds::Vector{Float64}`: Upper and Lower bounds for the reactive supply
+        at each bus.
 - `bus_type::Vector{PSY.ACBusTypes}`:
-        vector containing type of buses present in the system, ordered 
+        vector containing type of buses present in the system, ordered
         according to "bus_lookup".
 - `bus_magnitude::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus magnitudes, ordered according to 
+        "(b, t)" matrix containing the bus magnitudes, ordered according to
         "bus_lookup". b: number of buses, t: number of time period.
 - `bus_angles::Matrix{Float64}`:
-        "(b, t)" matrix containing the bus angles, ordered according to 
+        "(b, t)" matrix containing the bus angles, ordered according to
         "bus_lookup". b: number of buses, t: number of time period.
 - `branch_flow_values::Matrix{Float64}`:
-        "(br, t)" matrix containing the power flows, ordered according to 
+        "(br, t)" matrix containing the power flows, ordered according to
         "branch_lookup". br: number of branches, t: number of time period.
 - `timestep_map::Dict{Int, S}`:
-        dictonary mapping the number of the time periods (corresponding to the 
-        column number of the previosly mentioned matrices) and their actual 
-        names.
+        dictonary mapping the number of the time periods (corresponding to the
+        column number of the previosly mentioned matrices) and their names.
 - `valid_ix::Vector{Int}`:
-        vector containing the indeces related to those buses that are not slack
-        ones.
+        vector containing the indeces of not slack buses
 - `power_network_matrix::M`:
-        matrix used for the evaluation of either the power flows or bus angles, 
+        matrix used for the evaluation of either the power flows or bus angles,
         depending on the method considered.
 - `aux_network_matrix::N`:
-        matrix used for the evaluation of either the power flows or bus angles, 
+        matrix used for the evaluation of either the power flows or bus angles,
         depending on the method considered.
-
+- `neighbors::Vector{Set{Int}}`: Vector with the sets of adjacent buses.
 """
-struct PowerFlowData{M <: PNM.PowerNetworkMatrix, N, S <: Union{String, Char}}
+struct PowerFlowData{M <: PNM.PowerNetworkMatrix, N, S <: Union{String, String}}
     bus_lookup::Dict{Int, Int}
     branch_lookup::Dict{String, Int}
     bus_activepower_injection::Matrix{Float64}
     bus_reactivepower_injection::Matrix{Float64}
     bus_activepower_withdrawals::Matrix{Float64}
     bus_reactivepower_withdrawals::Matrix{Float64}
+    bus_reactivepower_bounds::Vector{Float64}
     bus_type::Vector{PSY.ACBusTypes}
     bus_magnitude::Matrix{Float64}
     bus_angles::Matrix{Float64}
@@ -63,10 +64,26 @@ struct PowerFlowData{M <: PNM.PowerNetworkMatrix, N, S <: Union{String, Char}}
     valid_ix::Vector{Int}
     power_network_matrix::M
     aux_network_matrix::N
+    neighbors::Vector{Set{Int}}
 end
 
 # AC Power Flow Data
 # TODO -> MULTI PERIOD: AC Power Flow Data
+function _calculate_neighbors(
+    Yb::PNM.Ybus{
+        Tuple{Vector{Int64}, Vector{Int64}},
+        Tuple{Dict{Int64, Int64}, Dict{Int64, Int64}},
+    },
+)
+    I, J, V = SparseArrays.findnz(Yb.data)
+    neighbors = [Set{Int}([i]) for i in 1:length(Yb.axes[1])]
+    for nz in eachindex(V)
+        push!(neighbors[I[nz]], J[nz])
+        push!(neighbors[J[nz]], I[nz])
+    end
+    return neighbors
+end
+
 """
 Function for the definition of the PowerFlowData strucure given the System
 data, number of time periods to consider and their names.
@@ -80,20 +97,23 @@ NOTE: use it for AC power flow computations.
         container storing the system data to consider in the PowerFlowData
         structure.
 - `timesteps::Int`:
-        number of time periods to consider in the PowerFlowData structure. It 
+        number of time periods to consider in the PowerFlowData structure. It
         defines the number of columns of the matrices used to store data.
         Default value = 1.
 - `timestep_names::Vector{String}`:
         names of the time periods defines by the argmunet "timesteps". Default
         value = String[].
+- `check_connectivity::Bool`:
+        Perform connectivity check on the network matrix. Default value = true.
 
 WARNING: functions for the evaluation of the multi-period AC OPF still to be implemented.
 """
 function PowerFlowData(
     ::ACPowerFlow,
-    sys::PSY.System,
+    sys::PSY.System;
     timesteps::Int = 1,
-    timestep_names::Vector{String} = String[])
+    timestep_names::Vector{String} = String[],
+    check_connectivity::Bool = true)
 
     # assign timestep_names
     # timestep names are then allocated in a dictionary to map matrix columns
@@ -106,7 +126,7 @@ function PowerFlowData(
     end
 
     # get data for calculations
-    power_network_matrix = PNM.Ybus(sys)
+    power_network_matrix = PNM.Ybus(sys; check_connectivity = check_connectivity)
 
     # get number of buses and branches
     n_buses = length(axes(power_network_matrix, 1))
@@ -119,14 +139,16 @@ function PowerFlowData(
     bus_lookup = power_network_matrix.lookup[2]
     branch_lookup =
         Dict{String, Int}(PSY.get_name(b) => ix for (ix, b) in enumerate(branches))
+
+    # TODO: bus_type might need to also be a Matrix since the type can change for a particular scenario
     bus_type = Vector{PSY.ACBusTypes}(undef, n_buses)
     bus_angles = zeros(Float64, n_buses)
     bus_magnitude = zeros(Float64, n_buses)
     temp_bus_map = Dict{Int, String}(
         PSY.get_number(b) => PSY.get_name(b) for b in PSY.get_components(PSY.Bus, sys)
     )
 
-    for (ix, bus_no) in bus_lookup
+    for (bus_no, ix) in bus_lookup
         bus_name = temp_bus_map[bus_no]
         bus = PSY.get_component(PSY.Bus, sys, bus_name)
         bus_type[ix] = PSY.get_bustype(bus)
@@ -185,14 +207,16 @@ function PowerFlowData(
         bus_reactivepower_injection_1,
         bus_activepower_withdrawals_1,
         bus_reactivepower_withdrawals_1,
+        Vector{Float64}(),
         bus_type,
         bus_magnitude_1,
         bus_angles_1,
         branch_flow_values_1,
-        Dict(zip([1], "1")),
+        Dict(1 => "1"),
         setdiff(1:n_buses, ref_bus_positions),
         power_network_matrix,
         nothing,
+        _calculate_neighbors(power_network_matrix),
     )
 end
 
@@ -205,24 +229,27 @@ NOTE: use it for DC power flow computations.
 
 # Arguments:
 - `::DCPowerFlow`:
-        use DCPowerFlow() to store the ABA matrix as power_network_matrix and 
+        use DCPowerFlow() to store the ABA matrix as power_network_matrix and
         the BA matrix as aux_network_matrix.
 - `sys::PSY.System`:
         container storing the system data to consider in the PowerFlowData
         structure.
 - `timesteps::Int`:
-        number of time periods to consider in the PowerFlowData structure. It 
+        number of time periods to consider in the PowerFlowData structure. It
         defines the number of columns of the matrices used to store data.
         Default value = 1.
 - `timestep_names::Vector{String}`:
         names of the time periods defines by the argmunet "timesteps". Default
         value = String[].
+- `check_connectivity::Bool`:
+        Perform connectivity check on the network matrix. Default value = true.
 """
 function PowerFlowData(
     ::DCPowerFlow,
-    sys::PSY.System,
+    sys::PSY.System;
     timesteps::Int = 1,
-    timestep_names::Vector{String} = String[])
+    timestep_names::Vector{String} = String[],
+    check_connectivity::Bool = true)
 
     # assign timestep_names
     # timestep names are then allocated in a dictionary to map matrix columns
@@ -293,8 +320,6 @@ function PowerFlowData(
     bus_angles_1 = deepcopy(init_1)
     branch_flow_values_1 = deepcopy(init_2)
 
-    # @show size(bus_activepower_injection_1)
-
     # initial values related to first timestep allocated in the first column
     bus_activepower_injection_1[:, 1] .= bus_activepower_injection
     bus_reactivepower_injection_1[:, 1] .= bus_reactivepower_injection
@@ -311,6 +336,7 @@ function PowerFlowData(
         bus_reactivepower_injection_1,
         bus_activepower_withdrawals_1,
         bus_reactivepower_withdrawals_1,
+        Vector{Float64}(),
         bus_type,
         bus_magnitude_1,
         bus_angles_1,
@@ -319,6 +345,7 @@ function PowerFlowData(
         setdiff(1:n_buses, aux_network_matrix.ref_bus_positions),
         power_network_matrix,
         aux_network_matrix,
+        Vector{Set{Int}}(),
     )
 end
 
@@ -331,13 +358,13 @@ NOTE: use it for DC power flow computations.
 
 # Arguments:
 - `::PTDFDCPowerFlow`:
-        use PTDFDCPowerFlow() to store the PTDF matrix as power_network_matrix 
+        use PTDFDCPowerFlow() to store the PTDF matrix as power_network_matrix
         and the ABA matrix as aux_network_matrix.
 - `sys::PSY.System`:
         container storing the system data to consider in the PowerFlowData
         structure.
 - `timesteps::Int`:
-        number of time periods to consider in the PowerFlowData structure. It 
+        number of time periods to consider in the PowerFlowData structure. It
         defines the number of columns of the matrices used to store data.
         Default value = 1.
 - `timestep_names::Vector{String}`:
@@ -346,9 +373,10 @@ NOTE: use it for DC power flow computations.
 """
 function PowerFlowData(
     ::PTDFDCPowerFlow,
-    sys::PSY.System,
+    sys::PSY.System;
     timesteps::Int = 1,
-    timestep_names::Vector{String} = String[])
+    timestep_names::Vector{String} = String[],
+    check_connectivity::Bool = true)
 
     # assign timestep_names
     # timestep names are then allocated in a dictionary to map matrix columns
@@ -437,6 +465,7 @@ function PowerFlowData(
         bus_reactivepower_injection_1,
         bus_activepower_withdrawals_1,
         bus_reactivepower_withdrawals_1,
+        Vector{Float64}(),
         bus_type,
         bus_magnitude_1,
         bus_angles_1,
@@ -445,6 +474,7 @@ function PowerFlowData(
         setdiff(1:n_buses, aux_network_matrix.ref_bus_positions),
         power_network_matrix,
         aux_network_matrix,
+        Vector{Set{Int}}(),
     )
 end
 
@@ -457,13 +487,13 @@ NOTE: use it for DC power flow computations.
 
 # Arguments:
 - `::PTDFDCPowerFlow`:
-        use vPTDFDCPowerFlow() to store the Virtual PTDF matrix as 
+        use vPTDFDCPowerFlow() to store the Virtual PTDF matrix as
         power_network_matrix and the ABA matrix as aux_network_matrix.
 - `sys::PSY.System`:
         container storing the system data to consider in the PowerFlowData
         structure.
 - `timesteps::Int`:
-        number of time periods to consider in the PowerFlowData structure. It 
+        number of time periods to consider in the PowerFlowData structure. It
         defines the number of columns of the matrices used to store data.
         Default value = 1.
 - `timestep_names::Vector{String}`:
@@ -472,9 +502,10 @@ NOTE: use it for DC power flow computations.
 """
 function PowerFlowData(
     ::vPTDFDCPowerFlow,
-    sys::PSY.System,
+    sys::PSY.System;
     timesteps::Int = 1,
-    timestep_names::Vector{String} = String[])
+    timestep_names::Vector{String} = String[],
+    check_connectivity::Bool = true)
 
     # assign timestep_names
     # timestep names are then allocated in a dictionary to map matrix columns
@@ -563,6 +594,7 @@ function PowerFlowData(
         bus_reactivepower_injection_1,
         bus_activepower_withdrawals_1,
         bus_reactivepower_withdrawals_1,
+        Vector{Float64}(),
         bus_type,
         bus_magnitude_1,
         bus_angles_1,
@@ -571,5 +603,6 @@ function PowerFlowData(
         setdiff(1:n_buses, aux_network_matrix.ref_bus_positions),
         power_network_matrix,
         aux_network_matrix,
+        Vector{Set{Int}}(),
     )
 end

src/PowerFlows.jl

@@ -1,7 +1,7 @@
 module PowerFlows
 
 export solve_powerflow
-export solve_powerflow!
+export solve_ac_powerflow!
 export PowerFlowData
 export DCPowerFlow
 export ACPowerFlow
@@ -14,7 +14,6 @@ import PowerSystems
 import LinearAlgebra
 import NLsolve
 import SparseArrays
-import SharedArrays: SharedArray
 import InfrastructureSystems
 import PowerNetworkMatrices
 import SparseArrays: SparseMatrixCSC
@@ -28,6 +27,8 @@ include("definitions.jl")
 include("powerflow_types.jl")
 include("PowerFlowData.jl")
 include("solve_dc_powerflow.jl")
+include("ac_power_flow.jl")
+include("ac_power_flow_jacobian.jl")
 include("nlsolve_ac_powerflow.jl")
 include("post_processing.jl")