Basic Data Utilities in PowerModels.jl---some more methods?

Specific Domains Optimization (Mathematical)
, ,
1

This post is a feature request related to Basic Data Utilities · PowerModels.

Although this sparse incidence matrix (denoted A) is useful, sometimes we may need another more succinct matrix—I call it bft, branch_from_to matrix—can be derived as follows

function A_2_bft(A) # from the incidence matrix to a B-by-2 branch_from_to matrix
    B, N = size(A)
    return bft = [findfirst(x -> x == (c == 1 ? -1 : 1), view(A, r, :)) for r in 1:B, c in 1:2]
end

Note that this matrix, aside from being more succinct, it is less ambiguous (-1, +1 in A is stipulated by the developer).
Could PowerModels.jl add a method that retrieve bft matrix directly?
And the conversion APIs A_2_bft, bft_2_A?
I would write (the +1, -1 might be swapped due to different conventions)

function bft_2_A(bft)
    B, N = size(bft, 1), maximum(bft)
    return SparseArrays.sparse([Vector(1:B); Vector(1:B)], vec(bft), [-ones(Int, B); ones(Int, B)], B, N)
end

We may further check the normality of a bft matrix

function assert_bft_is_normal(bft)
    minimum(bft) != 1 && error("smallest node label isn't 1")
    isempty(setdiff(1:maximum(bft), vec(bft))) || error("there is a gap between node labels")
end
2

@WalterMadelim, thanks for you interest in PowerModels!

First, it is best to make feature requests in the PowerModels github issues so we can debate and track them all in a common place.

Regarding this feature specifically. It am not immediately convinced this is an essential feature for the basic data utilities in PowerModels. The objective of basic data utilities was to provide a small collection of features that would support teaching of the most common power system analysis workflows. I worked closely with a few professors in the US power system community to determine what should be included or not. The incidence matrix as defined was determined as an important standard to support, as it supports workflows for computing PTDF and LODF matrices). I can imagine that your proposed branch_from_to matrix is useful, but if it is not a widely adopted standing in the power system classroom, then it is not a good fit for the basic data utilities of PowerModels.

1 Like
3

I see.

Are you aware that if there is a tool that can draw a node-branch diagram of the power system, given a specific incidence matrix A, or equivalently the bft?
I’ve just tried AI, the diagrams they generate are dishonest and unusable.

4

Possibly GitHub - WISPO-POP/PowerPlots.jl: Functions plot PowerModels networks is what you are looking for?

5

Not exactly. It totally depends on the established data format in those “cases(.m)”.
It would be better if it only takes an incidence matrix as input and draw the associated diagram.
Since I may make some custom modifications which not necessarily follow the data format of those "case.m"s.

"" 
# import LinearAlgebra.dot as dot
# import JuMP, Gurobi
import SparseArrays, PowerModels # ✅ branch numbering is more immediate, and node number lags behind

function relabel!(A, o, n) return @. A[A == o] = n end;
function A_2_bft(A) # from the incidence matrix to a B-by-2 branch_from_to matrix
    B, N = size(A)
    return bft = [findfirst(x -> x == (c == 1 ? -1 : 1), view(A, r, :)) for r in 1:B, c in 1:2]
end;
function bft_2_A(bft)
    B, N = size(bft, 1), maximum(bft)
    return SparseArrays.sparse([Vector(1:B); Vector(1:B)], vec(bft), [-ones(Int, B); ones(Int, B)], B, N)
end;
function assert_A_has_no_zero_col(A)
    for c in eachcol(A)
        isempty(SparseArrays.findnz(c)[1]) && error("A has a zero col")
    end
end;
function get_single_node_vec(A) return [n for n in 1:size(A, 2) if sum(abs.(view(A, :, n))) == 1] end;
function assert_is_A_row_normal(A)
    B, N = size(A)
    for b in 1:B # Make sure that each branch goes from one node to another different node.
        ns, vs = SparseArrays.findnz(A[b, :])
        vs == [1, -1] && continue
        vs == [-1, 1] && continue
        error("A is row-abnormal")
    end
end;
function assert_is_A_row_normal!(A)
    B, N = size(A)
    for b in 1:B
        ns, vs = SparseArrays.findnz(A[b, :])
        if vs[1] == 1
            A[b, :] *= -1
        end
    end
end;
function assert_is_A_row_normal2(A)
    B, N = size(A)
    for b in 1:B # Make sure that each branch goes from one node to another different node.
        ns, vs = SparseArrays.findnz(A[b, :])
        vs == [-1, 1] && continue
        error("2: A is row-abnormal")
    end
end;
function rectify_ref_dir_of_branches!(A)
    assert_is_A_row_normal(A)
    assert_is_A_row_normal!(A)
    assert_is_A_row_normal2(A)
end;
function get_b_on_of_single_node(bft, n)
    b, c = findfirst(x -> x == n, bft).I # one and only one
    on = bft[b, 3 - c]
    return b, on # branch, the other node
end;
function get_parallel_branches(bft)
    B, pbs = size(bft, 1), zeros(Int, 1, 2) # each row is a parallel branch pair
    for b in 1:B, ♭ in b+1:B
        (bft[♭, 1] == bft[b, 1] && bft[♭, 2] == bft[b, 2]) && (pbs = [pbs; b ♭])
    end
    size(pbs, 1) > 1 && return pbs[2:end, :]
    error("there is no parallel branches found")
end;
function get_pure_joint_bus(A, nvg, nvl)
    v = Int[]
    for n in setdiff(1:size(A, 2), union(nvg, nvl)) # loop over those "bare" nodes
        d = sum(abs.(view(A, :, n)))
        if d >= 3
            continue
        elseif d == 2
            push!(v, n)
        else
            error("node degree ≤ 1")
        end
    end
    return v
end;
function assert_is_connected(A)
    bft = A_2_bft(A);
    B, N = size(A)
    pb = Vector(1:B) # primal bs
    sn = [1] # subnet
    for ite in 1:B+1
        progress = false
        for (i, b) in enumerate(pb) # take out branch b
            if bft[b, 1] in sn
                on = bft[b, 2] # the other node
                on ∉ sn && push!(sn, on)
            elseif bft[b, 2] in sn
                on = bft[b, 1] # the other node
                on ∉ sn && push!(sn, on)
            else
                continue
            end
            popat!(pb, i); progress = true; break # fathom branch b
        end
        if progress == false
            error("ite = $ite. All the rest branches are not connected to the subnet being investigated. Check the current subnet")
        end
        1:N ⊆ sn && return # The graph is proved to be connected
    end
    error("here shouldn't be reached")
end;

# initial raw data
PowerModels.silence(); ⅁ = PowerModels.make_basic_network(PowerModels.parse_file("data/case118.m"));
nvg = [⅁["gen"]["$g"]["gen_bus"] for g in 1:length(⅁["gen"])]; # 🟠 a bus_index vector that has generators
nvl = [⅁["load"]["$l"]["load_bus"] for l in 1:length(⅁["load"])]; # 🟠 a bus index vector that has loads
A = -PowerModels.calc_basic_incidence_matrix(⅁);
𝐑 = [⅁["branch"]["$b"]["rate_a"] for b in 1:size(A, 1)]; # 🟠

assert_A_has_no_zero_col(A)
rectify_ref_dir_of_branches!(A)
single_node_vec = get_single_node_vec(A)

@assert issubset(single_node_vec, union(nvg, nvl))
issubset(single_node_vec, intersect(nvg, nvl)) || @info "the current topology can be enhanced"
# analyze the local structure
real_single_node_vec = setdiff(single_node_vec, intersect(nvl, nvg)) # these are single nodes to be removed, as they are simple
load_single_node_vec, gen_single_node_vec = intersect(real_single_node_vec, nvl), intersect(real_single_node_vec, nvg)
bft = A_2_bft(A); # 🟠
f = n -> get_b_on_of_single_node(bft, n)[1];
b_vec = f.(real_single_node_vec); # these branches are to be deleted!
f = n -> get_b_on_of_single_node(bft, n)[2];
load_single_on_vec, gen_single_on_vec = f.(load_single_node_vec), f.(gen_single_node_vec); # these are the end nodes after the action! to be done
# actions here!
nvl = union(setdiff(nvl, load_single_node_vec), load_single_on_vec);
nvg = union(setdiff(nvg, gen_single_node_vec), gen_single_on_vec);
ø = setdiff(1:size(bft, 1), b_vec); bft, 𝐑 = bft[ø, :], 𝐑[ø, :];
@info "these node_labels are absent: $(setdiff(1:maximum(bft), bft)), current max node_label = $(maximum(bft))"
let # 3 br-node deleted!
    relabel!(nvl, 115, 10);
    relabel!(nvg, 115, 10);
    relabel!(bft, 115, 10);
    relabel!(nvl, 116, 87);
    relabel!(nvg, 116, 87);
    relabel!(bft, 116, 87);
    relabel!(nvl, 118, 111);
    relabel!(nvg, 118, 111);
    relabel!(bft, 118, 111);
end;
A = bft_2_A(bft);
assert_A_has_no_zero_col(A)
rectify_ref_dir_of_branches!(A)
single_node_vec = get_single_node_vec(A)
@assert issubset(single_node_vec, union(nvg, nvl))
issubset(single_node_vec, intersect(nvg, nvl)) || @info "the current topology can be enhanced"
# analyze the local structure
real_single_node_vec = setdiff(single_node_vec, intersect(nvl, nvg)) # these are single nodes to be removed, as they are simple
load_single_node_vec, gen_single_node_vec = intersect(real_single_node_vec, nvl), intersect(real_single_node_vec, nvg)
bft = A_2_bft(A);
f = n -> get_b_on_of_single_node(bft, n)[1];
b_vec = f.(real_single_node_vec); # these branches are to be deleted!
f = n -> get_b_on_of_single_node(bft, n)[2];
load_single_on_vec, gen_single_on_vec = f.(load_single_node_vec), f.(gen_single_node_vec); # these are the end nodes after the action! to be done
# actions here!
nvg = union(setdiff(nvg, gen_single_node_vec), gen_single_on_vec);
ø = setdiff(1:size(bft, 1), b_vec); bft, 𝐑 = bft[ø, :], 𝐑[ø, :];
@info "these node_labels are absent: $(setdiff(1:maximum(bft), bft)), current max node_label = $(maximum(bft))"
let
    relabel!(nvl, 114, 9);
    relabel!(nvg, 114, 9);
    relabel!(bft, 114, 9);
end;
A = bft_2_A(bft);
assert_A_has_no_zero_col(A)
rectify_ref_dir_of_branches!(A)
single_node_vec = get_single_node_vec(A)
@assert issubset(single_node_vec, union(nvg, nvl))
issubset(single_node_vec, intersect(nvg, nvl)) || @info "the current topology can be enhanced"
bft = A_2_bft(A); # ✅ until here, we've finished deleting real_single_nodes
pbs = get_parallel_branches(bft);
for r in 1:size(pbs, 1)
    i, j = view(pbs, r, :)
    𝐑[i] += 𝐑[j]
end;
rsl = setdiff(1:size(bft, 1), view(pbs, :, 2)); # remnant rows
𝐑, bft = 𝐑[rsl], bft[rsl, :]; # remove! the redundant parallel branches
A = bft_2_A(bft);

nvi = get_pure_joint_bus(A, nvg, nvl); # a bus_index vector that is isolated (to be eliminated)
let
    @info "Visualize the nodes to be merged"
    for n in nvi
        (b1, c1), (b2, c2) = getproperty.(findall(x -> x == n, bft), :I)
        n1c = bft[b1, 3 - c1]
        n2c = bft[b2, 3 - c2]
        println("$n1c --($b1)-- $n --($b2)-- $n2c")
    end
end;
# add! new branches and delete! old branches
bft = [bft; 59 64]; push!(𝐑, min(𝐑[88], 𝐑[89]));
bft = [bft; 68 80]; push!(𝐑, min(𝐑[119], 𝐑[120]));
old_br_ind_vec = [88, 89, 119, 120];
deleteat!(𝐑, old_br_ind_vec);
bft = bft[setdiff(1:size(bft, 1), old_br_ind_vec), :];
let
    relabel!(bft, 112, 63)
    relabel!(bft, 113, 81)
    relabel!(nvl, 112, 63)
    relabel!(nvl, 113, 81)
    relabel!(nvg, 112, 63)
    relabel!(nvg, 113, 81)
end;
A = bft_2_A(bft);
assert_A_has_no_zero_col(A)
rectify_ref_dir_of_branches!(A)
bft = A_2_bft(A);
assert_is_connected(A);
sort!(nvg); sort!(nvl);
@info "Till here, we have derived a basic network that has none of 1. elongating line 2. parallel line 3. pure joint bus"
6

I do not know the Julia graph libraries well. @odow, do you have a suggestion on this topic?

7

Thanks for your willing to help. I’ve figured out a way.
This is the code to generate bft describing the graph.

import JuMP, Gurobi
import LinearAlgebra.dot as dot
import SparseArrays.sparse as sparse
import Random
function adjnode(n)
    nv = Int[]
    for ci in findall(x -> x == n, bft)
        b, c = ci.I
        on = bft[b, 3 - c] # the other node on branch b
        if on in nv
            error("in Function adjnode(), shouldn't happen")
        else
            push!(nv, on)
        end
    end
    return sort(nv)
end
function search_extreme_cases_code() # 📘 invoke this code on demand
    sd_vec, obj_vec = Int[], Float64[]
    while true
        sd = rand(Int) # seed
        Random.seed!(sd)
        JuMP.delete(ø, c)
        LD = rand(0.03:7e-15:1) * lrand(L) .* R_lnv; # randomly generate a load vector at an instance
        c = JuMP.@constraint(ø, [n in 1:N], ı_at(n) + dot( 𝐟 , view(A, :, n)) == w_at(n)); # 📘 KCL
        JuMP.optimize!(ø)
        if JuMP.termination_status(ø) == JuMP.OPTIMAL
            push!(sd_vec, sd)
            push!(obj_vec, JuMP.objective_value(ø))
            l = length(sd_vec)
            print("\rl = $l")
            l == 1 && break
        end
    end
    while true
        sd = rand(Int) # seed
        Random.seed!(sd)
        JuMP.delete(ø, c)
        LD = rand(0.03:7e-15:1) * lrand(L) .* R_lnv; # randomly generate a load vector at an instance
        c = JuMP.@constraint(ø, [n in 1:N], ı_at(n) + dot( 𝐟 , view(A, :, n)) == w_at(n)); # 📘 KCL
        JuMP.optimize!(ø)
        if JuMP.termination_status(ø) == JuMP.OPTIMAL
            o = JuMP.objective_value(ø)
            if o < obj_vec[1]
                pushfirst!(sd_vec, sd)
                pushfirst!(obj_vec, o)
            elseif o > obj_vec[end]
                push!(sd_vec, sd)
                push!(obj_vec, o)
            else
                continue 
            end
            l = length(sd_vec)
            print("\rl = $l")
            l == 20 && break
        end
    end
end
function lrand() return sin(rand(0.1:7e-15:1.57)) end;
function lrand(N) return sin.(rand(0.1:7e-15:1.57, N)) end;
function bft_2_A(bft)
    B, N = size(bft, 1), maximum(bft)
    return sparse([Vector(1:B); Vector(1:B)], vec(bft), [-ones(Int, B); ones(Int, B)], B, N)
end;
function load_111_data()
    R = [10.76, 25.35, 137.5, 10.15, 20.32, 52.75, 42.08, 15.63, 15.79, 54.85, 17.45, 6.72, 32.03, 14.7, 15.21, 4.4, 5.51, 13.06, 24.61, 6.05, 21.62, 22.23, 9.39, 27.28, 12.94, 11.32, 6.91, 22.02, 13.79, 29.41, 6.77, 12.82, 11.56, 28.96, 22.21, 13.01, 6.88, 32.27, 9.4, 10.92, 14.24, 8.64, 4.35, 107.59, 22.07, 7.59, 39.87, 115.33, 29.96, 10.15, 6.31, 20.73, 17.77, 22.11, 5.88, 7.96, 4.44, 6.49, 12.1, 7.95, 8.48, 5.67, 17.19, 6.79, 5.67, 20.97, 14.08, 7.73, 18.06, 6.67, 9.0, 7.47, 15.46, 113.06, 70.8, 10.96, 7.91, 10.96, 14.73, 4.79, 8.71, 5.09, 7.57, 7.32, 81.68, 19.56, 29.19, 41.93, 11.35, 37.06, 24.0, 5.03, 9.38, 30.37, 10.81, 69.97, 3.87, 3.32, 30.37, 8.61, 2.73, 30.72, 5.56, 6.06, 24.31, 8.13, 7.63, 8.74, 26.49, 7.27, 10.64, 5.38, 86.71, 44.94, 32.17, 15.58, 12.44, 29.32, 7.69, 7.29, 15.86, 8.79, 10.81, 6.43, 15.49, 16.72, 12.86, 28.73, 8.45, 12.68, 6.8, 14.68, 24.77, 6.06, 20.27, 12.35, 11.81, 10.16, 5.33, 3.72, 18.52, 19.61, 12.46, 6.13, 13.49, 8.7, 19.63, 9.8, 20.47, 5.38, 6.81, 6.57, 4.74, 28.75, 19.9, 5.9, 14.98, 5.9, 36.65, 6.06, 13.85, 16.38, 35.72, 5.3, 17.93, 14.81, 105.49, 276.46, 22.37, 19.78, 29.11, 30.37]
    gnv = [1, 4, 6, 8, 12, 15, 18, 19, 24, 25, 26, 27, 31, 32, 34, 36, 40, 42, 46, 49, 54, 55, 56, 59, 61, 62, 63, 65, 66, 69, 70, 72, 73, 74, 76, 77, 80, 81, 85, 86, 87, 89, 90, 91, 92, 99, 100, 103, 104, 105, 107, 110]
    lnv = [1, 2, 3, 4, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 27, 28, 29, 31, 32, 33, 34, 35, 36, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 62, 63, 66, 67, 70, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111]
    bft = [1 2; 1 3; 4 5; 3 5; 5 6; 6 7; 5 8; 4 11; 5 11; 11 12; 2 12; 3 12; 7 12; 11 13; 12 14; 13 15; 14 15; 12 16; 15 17; 16 17; 17 18; 18 19; 19 20; 15 19; 20 21; 21 22; 22 23; 23 24; 23 25; 25 26; 25 27; 27 28; 28 29; 17 30; 8 30; 26 30; 17 31; 29 31; 23 32; 31 32; 27 32; 15 33; 19 34; 35 36; 35 37; 33 37; 34 36; 34 37; 37 38; 37 39; 37 40; 30 38; 39 40; 40 41; 40 42; 41 42; 43 44; 34 43; 44 45; 45 46; 46 47; 46 48; 47 49; 42 49; 45 49; 48 49; 49 50; 49 51; 51 52; 52 53; 53 54; 49 54; 54 55; 54 56; 55 56; 56 57; 50 57; 56 58; 51 58; 54 59; 56 59; 55 59; 59 60; 59 61; 60 61; 60 62; 61 62; 61 64; 38 65; 64 65; 49 66; 62 66; 62 67; 65 66; 66 67; 65 68; 47 69; 49 69; 68 69; 69 70; 24 70; 70 71; 24 72; 71 72; 71 73; 70 74; 70 75; 69 75; 74 75; 76 77; 69 77; 75 77; 77 78; 78 79; 77 80; 79 80; 77 82; 82 83; 83 84; 83 85; 84 85; 85 86; 85 88; 85 89; 88 89; 89 90; 90 91; 89 92; 91 92; 92 93; 92 94; 93 94; 94 95; 80 96; 82 96; 94 96; 80 97; 80 98; 80 99; 92 100; 94 100; 95 96; 96 97; 98 100; 99 100; 100 101; 92 102; 101 102; 100 103; 100 104; 103 104; 103 105; 100 106; 104 105; 105 106; 105 107; 105 108; 106 107; 108 109; 103 110; 109 110; 63 110; 17 81; 32 81; 9 32; 10 27; 9 10; 68 87; 75 111; 76 111; 59 64; 68 80]
    return R, gnv, lnv, bft # RateA, generator_node_vector, load_node_vector, branch_from_to
end;
function ı_at(n) # return the power injection Given a node index
    i = findfirst(x -> x == n, gnv) # injection index, NOT a generator index
    b = isnothing(i)
    if b
        return !b # no injection
    else
        return ı[i] # a VariableRef
    end
end;
function w_at(n) # return the power withdrawal Given a node index
    l = findfirst(x -> x == n, lnv) # load index
    b = isnothing(l)
    if b
        return !b # no injection
    else
        return LD[l] # load data
    end
end;
R, gnv, lnv, bft = load_111_data(); A = bft_2_A(bft);

And this is the code that generate a graph that can be understood by human.


using Graphs, GraphPlot, CairoMakie
function clockwise_rotate(A, t)
    c, s = cos(t), sin(t)
    [c s; -s c] * A
end
function relocate_as_mean!(xv, yv, n, n1, n2)
    xv[n] = (xv[n1] + xv[n2])/2
    yv[n] = (yv[n1] + yv[n2])/2
end
g = SimpleGraph(maximum(bft)); for (n, m) in eachrow(bft) # graph construction
    add_edge!(g, n, m)
end
function custom_modification!(xv, yv)
    relocate_as_mean!(xv, yv, 48, 66, 68)
    relocate_as_mean!(xv, yv, 47, 46, 69)
    relocate_as_mean!(xv, yv, 50, 49, 57)
    relocate_as_mean!(xv, yv, 92, 91, 94)
    relocate_as_mean!(xv, yv, 25, 23, 27)
    relocate_as_mean!(xv, yv, 9, 29, 32)
    relocate_as_mean!(xv, yv, 86, 83, 89)
    relocate_as_mean!(xv, yv, 79, 77, 97)
    relocate_as_mean!(xv, yv, 87, 69, 80)
    relocate_as_mean!(xv, yv, 51, 42, 50)
    relocate_as_mean!(xv, yv, 73, 71, 74)
    relocate_as_mean!(xv, yv, 43, 15, 44)
    relocate_as_mean!(xv, yv, 1, 2, 14)
    relocate_as_mean!(xv, yv, 108, 55, 107)
    relocate_as_mean!(xv, yv, 109, 108, 110)
    relocate_as_mean!(xv, yv, 63, 91, 110)
    relocate_as_mean!(xv, yv, 88, 86, 89)
    relocate_as_mean!(xv, yv, 90, 89, 91)
    relocate_as_mean!(xv, yv, 72, 23, 69)
end; Random.seed!(2527829590547944688);
xv, yv = spring_layout(g); # supported by the layout algorithm, crude
angle = 39; # choose an apt angle
temp = clockwise_rotate([xv'; yv'], π/180 * angle);
xv, yv = temp[1, :], temp[2, :];
custom_modification!(xv, yv);

line_colors = [
    :chocolate, # Yellow
    :gold1,     # Yellow
    :blue,      # blue 
    :cyan,      # blue
    :lime,      # green lemon
    :darkgreen, # green
    :olive,     # green
    :deeppink,  # Red
    :red4       # Red
];
label_size = 8

f = Figure(); ax = Axis(f[1, 1]); # 🖼️
Random.seed!(2)
for (n, m) in eachrow(bft) # Draw lines
    lines!(ax, [xv[n], xv[m]], [yv[n], yv[m]]; linewidth = 1, color = rand(line_colors))
end
for n in 1:maximum(bft) # Draw Points
    if n in gnv
        if n in lnv
            text!(ax, [xv[n]], [yv[n]]; color = :purple, text = "$n", fontsize = label_size)
        else
            text!(ax, [xv[n]], [yv[n]]; color = :red, text = "$n", fontsize = label_size)
        end
    elseif n in lnv
        text!(ax, [xv[n]], [yv[n]]; color = :skyblue, text = "$n", fontsize = label_size)
    else
        text!(ax, [xv[n]], [yv[n]]; color = :black, text = "$n", fontsize = label_size)
    end
end
hidedecorations!(ax)
f # After typing this in REPL, it will generate the picture

The picture derived is

image
image1916×1412 132 KB

The good ideas include

  1. use diverse line colors (I use 2 Yellow + 2 Blue + 2 Red + 3 Green)
  2. as a user, you need to do some custom relocations inevitably. First, Givens rotation can be employed to make the layout look upright. Second, adjust specific point locations, as I had done in custom_modification!.
1 Like
8

I don’t have any specific suggestions

9

Shameless plug: you might want to check out GraphMakie.jl, which is doing something quite similar to what you did here: take a graph and style it according to some rules, including node size, edge thickness, labels, colors, …

For bus placement its especially helpful to look at the pin keyword arguments of Stress and Spring layouts, which can be used to place some of the vertices manually, while the rest are layouted automaticially around and between the fixed nodes.

1 Like
10

Perhaps you know about it already, but sometimes for highly detailed graphs with not too many nodes I prefer GraphViz (https://graphviz.org/). It’s easy to build the DOT language string in Julia and send it to the command line tool to render it into SVG, PDF, etc. I did this for the multi-commodity network problem solved with JuMP here JuMP-ing with AlgebraicJulia II: A practical optimization model – AlgebraicJulia blog (see end of post for the visualization)

2 Likes