Incidence, BA and ABA matrices

In this tutorial the IncidenceMatrix, BA_matrix and ABA_matrix are presented. The methods used for their evaluation, as well as how data is stored is shown in the following subsections.

The matrices here presented are the building blocks for the computation of the PTDF and LODF matrices.

IncidenceMatrix

The PowerNetworkMatrices package defines the structure IncidenceMatrix, which store the Incidence Matrix of the considered system as well as the most relevant network data.

At first, the System data is loaded.

julia> using PowerNetworkMatrices
julia> using PowerSystemCaseBuilder
julia> import PowerNetworkMatrices as PNM
julia> import PowerSystemCaseBuilder as PSB
julia> sys = PSB.build_system(PSB.PSITestSystems, "c_sys5");[ Info: Loaded time series from storage file existing=/home/runner/.julia/packages/PowerSystemCaseBuilder/zW01F/data/serialized_system/8584b9e729c8aa68ee5405660c6258cde1f67ed3b68f114823707912e9a0d16c/c_sys5_time_series_storage.h5 new=/tmp/jl_pVhzFU compression=InfrastructureSystems.CompressionSettings(false, InfrastructureSystems.CompressionTypesModule.CompressionTypes.DEFLATE = 1, 3, true)

Then the Incidence Matrix is computed as follows:

julia> incidence_matrix = PNM.IncidenceMatrix(sys)[ Info: Finding subnetworks via iterative union find
IncidenceMatrix{Tuple{Vector{Tuple{Int64, Int64}}, Vector{Int64}}, Tuple{Dict{Tuple{Int64, Int64}, Int64}, Dict{Int64, Int64}}}
    Dimension 1, [(4, 5), (1, 2), (3, 4), (1, 4), (2, 3), (1, 5)]
    Dimension 2, [1, 2, 3, 4, 5]
And data with size (6, 5):
 ⋅   ⋅   ⋅   1  -1
 1  -1   ⋅   ⋅   ⋅
 ⋅   ⋅   1  -1   ⋅
 1   ⋅   ⋅  -1   ⋅
 ⋅   1  -1   ⋅   ⋅
 1   ⋅   ⋅   ⋅  -1

The incidence_matrix variable is a structure of type IncidenceMatrix. Getter functions are available to access additional information:

julia> # axis names: row and column names.
       # row names: tuples of the arcs (from bus number, to bus number)
       # column names: names of the buses
       PNM.get_axes(incidence_matrix)([(4, 5), (1, 2), (3, 4), (1, 4), (2, 3), (1, 5)], [1, 2, 3, 4, 5])
julia> # data: the incidence matrix data PNM.get_data(incidence_matrix)6×5 SparseArrays.SparseMatrixCSC{Int8, Int64} with 12 stored entries: ⋅ ⋅ ⋅ 1 -1 1 -1 ⋅ ⋅ ⋅ ⋅ ⋅ 1 -1 ⋅ 1 ⋅ ⋅ -1 ⋅ ⋅ 1 -1 ⋅ ⋅ 1 ⋅ ⋅ ⋅ -1
julia> # lookup: dictionary linking the arc tuples and bus numbers with the row # and column numbers, respectively. PNM.get_lookup(incidence_matrix)(Dict((4, 5) => 1, (1, 2) => 2, (3, 4) => 3, (1, 4) => 4, (2, 3) => 5, (1, 5) => 6), Dict(5 => 5, 4 => 4, 2 => 2, 3 => 3, 1 => 1))
julia> # ref_bus_positions: set containing the positions of the reference buses. # this represents the positions where to add the column of zeros. Please refer to the # example in the BA matrix for more details. PNM.get_ref_bus_position(incidence_matrix)1-element Vector{Int64}: 4

Note that the number of columns is lower than the actual number of system buses since the column related to the reference bus is discarded.

BA_Matrix

The BA_Matrix is a structure containing the matrix coming from the product of the IncidenceMatrix and the diagonal matrix containing the impedance of the system's branches ("B" matrix).

The BA_Matrix is computed as follows:

julia> ba_matrix = PNM.BA_Matrix(sys)[ Info: Finding subnetworks via iterative union find
DC_BA_Matrix
    Dimension 1, [(4, 5), (1, 2), (3, 4), (1, 4), (2, 3), (1, 5)]
    Dimension 2, [1, 2, 3, 4, 5]
Note!! The data shown below corresponds to the transposed matrix.
And data with size (5, 6):
    ⋅     35.5872     ⋅     32.8947     ⋅       156.25
    ⋅    -35.5872     ⋅       ⋅       92.5926      ⋅
    ⋅       ⋅       33.67     ⋅      -92.5926      ⋅
  33.67     ⋅      -33.67  -32.8947     ⋅          ⋅
 -33.67     ⋅         ⋅       ⋅         ⋅      -156.25

Note that the axes order matches the IncidenceMatrix (arcs x buses), but for computational considerations the raw data is the transposed matrix.

The matrix data can similarly be accessed with getter functions:

julia> PNM.get_data(ba_matrix)5×6 SparseArrays.SparseMatrixCSC{Float64, Int64} with 12 stored entries:
    ⋅     35.5872     ⋅     32.8947     ⋅       156.25
    ⋅    -35.5872     ⋅       ⋅       92.5926      ⋅
    ⋅       ⋅       33.67     ⋅      -92.5926      ⋅
  33.67     ⋅      -33.67  -32.8947     ⋅          ⋅
 -33.67     ⋅         ⋅       ⋅         ⋅      -156.25

Note that the number of columns is lower than the actual number of system buses since the column related to the reference bus is discarded.

To add the column of zeros related to the reference bus, it is necessary to use the information from get_ref_bus_position.

julia> # assumes a single reference bus
       ref_bus_position = first(PNM.get_ref_bus_position(ba_matrix))4
julia> new_ba_matrix = hcat( ba_matrix.data[:, 1:(ref_bus_position - 1)], zeros(size(ba_matrix, 1), 1), ba_matrix.data[:, ref_bus_position:end], )5×7 SparseArrays.SparseMatrixCSC{Float64, Int64} with 12 stored entries: ⋅ 35.5872 ⋅ ⋅ 32.8947 ⋅ 156.25 ⋅ -35.5872 ⋅ ⋅ ⋅ 92.5926 ⋅ ⋅ ⋅ 33.67 ⋅ ⋅ -92.5926 ⋅ 33.67 ⋅ -33.67 ⋅ -32.8947 ⋅ ⋅ -33.67 ⋅ ⋅ ⋅ ⋅ ⋅ -156.25

However, trying to change the data field with a matrix of different dimension will result in an error.

julia> ba_matrix.data = new_ba_matrixERROR: setfield!: immutable struct of type BA_Matrix cannot be changed

ABA_Matrix

The ABA_Matrix is a structure containing the matrix coming from the product of the IncidenceMatrix and the BA_Matrix. It features the same fields as the IncidenceMatrix and the BA_Matrix, plus the K one. The field ABA_Matrix.K stores the LU factorization matrices (using the methods contained in the package KLU).

To evaluate the ABA_Matrix, the following command is sufficient:

julia> aba_matrix = ABA_Matrix(sys);[ Info: Finding subnetworks via iterative union find

By default the LU factorization matrices are not computed, leaving the K field empty:

julia> isnothing(aba_matrix.K)true

In case these are wanted, the keyword factorize must be true.

julia> aba_matrix_with_LU = ABA_Matrix(sys; factorize = true);[ Info: Finding subnetworks via iterative union find
julia> aba_matrix_with_LU.KKLU.KLUFactorization{Float64, Int64} L factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 7 stored entries: 1.0 ⋅ ⋅ ⋅ -0.695273 1.0 ⋅ ⋅ ⋅ ⋅ 1.0 ⋅ ⋅ -0.648697 -0.722365 1.0 U factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 7 stored entries: 1.0 -0.822715 ⋅ ⋅ ⋅ 0.427989 ⋅ -0.158354 ⋅ ⋅ 1.0 -0.733333 ⋅ ⋅ ⋅ 0.367542 F factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 0 stored entries: ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅

If the ABA_Matrix is already computed but the LU factorization was not performed, this can be done by considering the following command:

julia> aba_matrix.K
julia> aba_matrix = factorize(aba_matrix);
julia> aba_matrix.KKLU.KLUFactorization{Float64, Int64} L factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 7 stored entries: 1.0 ⋅ ⋅ ⋅ -0.695273 1.0 ⋅ ⋅ ⋅ ⋅ 1.0 ⋅ ⋅ -0.648697 -0.722365 1.0 U factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 7 stored entries: 1.0 -0.822715 ⋅ ⋅ ⋅ 0.427989 ⋅ -0.158354 ⋅ ⋅ 1.0 -0.733333 ⋅ ⋅ ⋅ 0.367542 F factor: 4×4 SparseArrays.SparseMatrixCSC{Float64, Int64} with 0 stored entries: ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅

The following command can then be used to check if the ABA_Matrix contains the LU factorization matrices:

julia> is_factorized(aba_matrix)true