(home)

Checking if an electric power system is radial with the admittance matrix

A power system is "radial" if there are no loops in the circuit; that is, it represents a tree graph. In the context of electric power systems, this means that the network admittance matrix YY has a certain structure that can be checked using its off-diagonal elements.

Simply put, a radial network has no electrical loops, which means that the number of non-zero upper and lower off-diagonal elements in the admittance matrix is limited. If the network admittance matrix is YCn×nY \in \mathbb{C}^{n \times n}, then the network is not radial if

ijLijn,where (i,j){(i,j):Lij0}, \sum_{i} \sum_j L_{ij} \geq n,\\ \text{where} \ (i,j) \in \{(i,j): L_{ij} \neq 0 \},

or alternatively,

ijUijn,where (i,j){(i,j):Uij0}. \sum_i \sum_j U_{ij} \geq n,\\ \text{where} \ (i,j) \in \{(i,j): U_{ij} \neq 0 \}.

This function tells you if a PowerModels.jl network model is radial or not.

using LinearAlgebra

"""
Given a network's admittance matrix Y, determine if that network is radial.
In other words, determine if there are no electrical loops.
    Params:
    Y: NxN Network admittance matrix
"""
function is_radial(Y::AbstractMatrix)
    
    n = size(Y)[1]

    #Upper and lower off-diagonal elements
    U(A::AbstractMatrix) = [A[i] for i in CartesianIndices(A) if i[1]>i[2]]
    L(A::AbstractMatrix) = [A[i] for i in CartesianIndices(A) if i[1]<i[2]]
    
    #Get the nonzero upper and lower off diagonal elements
    nz_upper = [1 for y_ij in U(Y) if y_ij != 0]
    nz_lower = [1 for y_ij in L(Y) if y_ij != 0]
    
    return !(sum(nz_upper)>n-1 || sum(nz_lower) >n-1)
end

#EXAMPLE: case4_dist admittance matrix
Y = [
    133.333-266.667im -66.6667+133.333im 0.0+0.0im -65.0407+130.081im;
    -66.6667+133.333im 133.333-266.667im -66.6667+133.333im 0.0+0.0im;
    0.0+0.0im -66.6667+133.333im  66.6667-133.333im  0.0+0.0im;
    -65.0407+130.081im 0.0+0.0im 0.0+0.0im 63.4543-126.909im
    ]

println(is_radial(Y))
true