# Code

This is a page for short code blocks that I personally think are nice.

More code is on my Github.

Work in progress! Check back in the future.

## Check if an electric power system is radial

A power system is "radial" if there are no loops in the circuit.

If we call $$L$$ and $$U$$ the lower and upper off-diagonal components the network admittance matrix $$Y \in \mathbb{C}^{n \times n}$$, then the network is not radial if

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

or alternatively,

$\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:
"""

n = size(Y)

#Upper and lower off-diagonal elements
U(A::AbstractMatrix) = [A[i] for i in CartesianIndices(A) if i>i]
L(A::AbstractMatrix) = [A[i] for i in CartesianIndices(A) if i<i]

#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

println(is_radial(Y))
true