diff --git a/LICENSE.md b/LICENSE.md index ac743ed..fa45b11 100644 --- a/LICENSE.md +++ b/LICENSE.md @@ -67,3 +67,26 @@ GraphPlot.jl >The above copyright notice and this permission notice shall be included in all copies or substantial portions of the >Software. > >THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE >WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR >COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR >OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +rust-sugiyama +>MIT License +> +>Copyright (c) 2024 paddison +> +>Permission is hereby granted, free of charge, to any person obtaining a copy +>of this software and associated documentation files (the "Software"), to deal +>in the Software without restriction, including without limitation the rights +>to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +>copies of the Software, and to permit persons to whom the Software is +>furnished to do so, subject to the following conditions: +> +>The above copyright notice and this permission notice shall be included in all +>copies or substantial portions of the Software. +> +>THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +>IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +>FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +>AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +>LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +>OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE +>SOFTWARE. diff --git a/docs/src/index.md b/docs/src/index.md index d7d80a0..b95b478 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -82,6 +82,30 @@ f, ax, p = graphplot(g, layout=layout) hidedecorations!(ax); hidespines!(ax); ax.aspect = DataAspect(); f #hide ``` +## Sugiyama Layered Layout +```@docs +Sugiyama +``` +### Example +```@example layouts +g = SimpleDiGraph(10) +for (s,d) in [(1,2),(2,3),(3,4),(3,5),(4,6),(4,7),(4,8),(4,9),(5,6),(5,7),(5,8),(5,9), + (6,10),(7,10),(8,10),(9,10)] + add_edge!(g, s, d) +end +layout = Sugiyama() +f, ax, p = graphplot(g, layout=layout) +hidedecorations!(ax); hidespines!(ax); ax.aspect = DataAspect(); f #hide +``` +Edges that span more than one rank are routed around the nodes in between +via internal "dummy" vertices; use `direction=:right` to flow left-to-right +instead of top-to-bottom: +```@example layouts +layout = Sugiyama(direction=:right) +f, ax, p = graphplot(g, layout=layout) +hidedecorations!(ax); hidespines!(ax); ax.aspect = DataAspect(); f #hide +``` + ## Spring/Repulsion Model ```@docs Spring diff --git a/src/NetworkLayout.jl b/src/NetworkLayout.jl index 257fb88..951ded1 100644 --- a/src/NetworkLayout.jl +++ b/src/NetworkLayout.jl @@ -239,5 +239,6 @@ include("spectral.jl") include("shell.jl") include("squaregrid.jl") include("align.jl") +include("sugiyama.jl") end diff --git a/src/sugiyama.jl b/src/sugiyama.jl new file mode 100644 index 0000000..1968343 --- /dev/null +++ b/src/sugiyama.jl @@ -0,0 +1,196 @@ +export Sugiyama, sugiyama + +include("sugiyama/graph.jl") +include("sugiyama/cycles.jl") +include("sugiyama/ranking.jl") +include("sugiyama/ordering.jl") +include("sugiyama/coordinates.jl") + +""" + Sugiyama(; kwargs...)(adj_matrix) + sugiyama(adj_matrix; kwargs...) + +Layered ("hierarchical") layout for directed graphs: vertices are grouped +into ranks, edge crossings between ranks are heuristically minimized, and +coordinates are assigned within each rank. Ranking and crossing +minimization follow Gansner, Koutsofios, North & Vo (1993, +[doi 10.1109/32.221135](https://doi.org/10.1109/32.221135)); coordinate +assignment follows Brandes & Köpf (2002, +[doi 10.1007/3-540-45848-4_3](https://doi.org/10.1007/3-540-45848-4_3)). + +Takes the adjacency matrix of a directed graph and returns coordinates of +the nodes. Cycles are broken by implicitly reversing edges; disconnected +components are laid out independently and placed side by side. + +## Keyword Arguments +- `Ptype=Float64`: Determines the output type `Point{2,Ptype}`. +- `nodesize=Float64[]`: Size of each node. Filled up with `ones` or + truncated to match the number of nodes. +- `nodespacing=1.0`: Minimum gap between neighboring nodes, both within a + rank and between ranks. +- `dummysize=0.0`: Width of the invisible "dummy" nodes used internally to + route edges spanning more than one rank. +- `minimum_length=1`: Minimum number of ranks every edge must span. +- `ranking_type=:networksimplex`: `:networksimplex` minimizes total edge + length (as in the papers above); `:longestpath`, `:up` and `:down` are + cheaper longest-path schedules (from both ends, from sources only, and + from sinks only, respectively). +- `crossing_minimization=:barycenter`: `:barycenter` or `:median` heuristic + used to reorder each rank. +- `transpose=true`: Follow up with a greedy pairwise-swap pass that further + reduces crossings, at the cost of runtime. +- `direction=:down`: Which way ranks flow: `:down`, `:up`, `:left` or + `:right`. + +This implementation is a Julia port of +[rust-sugiyama](https://github.com/paddison/rust-sugiyama). +""" +@addcall struct Sugiyama{Ptype,T} <: AbstractLayout{2,Ptype} + nodesize::Vector{T} + nodespacing::Float64 + dummysize::Float64 + minimum_length::Int + ranking_type::Symbol + crossing_minimization::Symbol + transpose::Bool + direction::Symbol +end + +function Sugiyama(; Ptype=Float64, + nodesize=Float64[], + nodespacing=1.0, + dummysize=0.0, + minimum_length=1, + ranking_type=:networksimplex, + crossing_minimization=:barycenter, + transpose=true, + direction=:down) + minimum_length >= 1 || throw(ArgumentError("minimum_length must be >= 1")) + direction in (:down, :up, :left, :right) || + throw(ArgumentError("direction must be one of :down, :up, :left, :right")) + ranking_type in (:networksimplex, :longestpath, :up, :down) || + throw(ArgumentError("ranking_type must be one of :networksimplex, :longestpath, :up, :down")) + crossing_minimization in (:barycenter, :median) || + throw(ArgumentError("crossing_minimization must be :barycenter or :median")) + Sugiyama{Ptype,eltype(nodesize)}(nodesize, Float64(nodespacing), Float64(dummysize), + Int(minimum_length), ranking_type, crossing_minimization, + transpose, direction) +end + +function layout(algo::Sugiyama{Ptype,T}, adj_matrix::AbstractMatrix) where {Ptype,T} + n = assertsquare(adj_matrix) + positions = Vector{Point{2,Ptype}}(undef, n) + n == 0 && return positions + + nodesize = ones(Float64, n) + for i in 1:min(n, length(algo.nodesize)) + nodesize[i] = Float64(algo.nodesize[i]) + end + + edgelist = Tuple{Int,Int}[] + for j in 1:n, i in 1:n + i != j && !iszero(adj_matrix[i, j]) && push!(edgelist, (i, j)) + end + + xoffset = 0.0 + for comp in _weakly_connected_components(n, edgelist) + m = length(comp) + localid = Dict(v => k for (k, v) in enumerate(comp)) + + sub = SugiGraph() + for v in comp + _add_vertex!(sub; width=nodesize[v] + algo.nodespacing, height=nodesize[v] + algo.nodespacing) + end + for (i, j) in edgelist + (haskey(localid, i) && haskey(localid, j)) || continue + _add_edge!(sub, localid[i], localid[j]) + end + + xs, ys = layout_component!(sub, algo) + + minx, maxx = extrema(view(xs, 1:m)) + for (k, v) in enumerate(comp) + positions[v] = to_point(Ptype, xs[k] - minx + xoffset, ys[k], algo.direction) + end + xoffset += (maxx - minx) + algo.nodespacing + end + return positions +end + +function to_point(::Type{Ptype}, x::Float64, y::Float64, direction::Symbol) where {Ptype} + if direction === :down + Point{2,Ptype}(x, -y) + elseif direction === :up + Point{2,Ptype}(x, y) + elseif direction === :right + Point{2,Ptype}(y, x) + else # :left + Point{2,Ptype}(-y, x) + end +end + +"""Weakly connected components of the graph given by `n` vertices `1:n` and +directed edge list `edgelist`, as sorted vectors of (global) vertex ids.""" +function _weakly_connected_components(n::Int, edgelist::Vector{Tuple{Int,Int}}) + adj = [Int[] for _ in 1:n] + for (i, j) in edgelist + push!(adj[i], j) + push!(adj[j], i) + end + visited = falses(n) + comps = Vector{Int}[] + for start in 1:n + visited[start] && continue + comp = Int[] + stack = [start] + visited[start] = true + while !isempty(stack) + v = pop!(stack) + push!(comp, v) + for w in adj[v] + if !visited[w] + visited[w] = true + push!(stack, w) + end + end + end + push!(comps, sort!(comp)) + end + return comps +end + +"""Run all four layout phases on a single weakly-connected component `g` +(vertices `1:m` are the "real" input vertices; the pipeline may append +dummy vertices after that). Returns `(xs, ys)` covering every vertex of the +(possibly grown) graph.""" +function layout_component!(g::SugiGraph, algo::Sugiyama) + remove_cycles!(g) + rank!(g, algo.minimum_length, algo.ranking_type) + insert_dummy_vertices!(g, algo.minimum_length, algo.dummysize) + layers = ordering(g, algo.crossing_minimization, algo.transpose) + + layouts = create_layouts(g, layers) + align_to_smallest_width_layout!(layouts) + xcoords = calculate_relative_coords(layouts) + + n = _nv(g) + xs = [xcoords[v] for v in 1:n] + xs .-= minimum(xs) + + rank_to_max_height = Dict{Int,Float64}() + for v in 1:n + r = g.verts[v].rank + rank_to_max_height[r] = max(get(rank_to_max_height, r, 0.0), g.verts[v].height) + end + ranks_sorted = sort!(collect(keys(rank_to_max_height))) + rank_to_y = Dict{Int,Float64}() + top = -rank_to_max_height[ranks_sorted[1]] * 0.5 + for r in ranks_sorted + mh = rank_to_max_height[r] + rank_to_y[r] = top + mh * 0.5 + top += mh + end + ys = [rank_to_y[g.verts[v].rank] for v in 1:n] + + return xs, ys +end diff --git a/src/sugiyama/coordinates.jl b/src/sugiyama/coordinates.jl new file mode 100644 index 0000000..ee3daff --- /dev/null +++ b/src/sugiyama/coordinates.jl @@ -0,0 +1,253 @@ +# Phase 3: assign x/y coordinates given the ranks and per-rank ordering, +# using the horizontal-alignment method of Brandes & Köpf (2002), "Fast and +# Simple Horizontal Coordinate Assignment". +# +# The idea: align each vertex with (at most) one "median" neighbor per rank +# to build vertical chains ("blocks"), then compact each block as far left +# as it can go without overlapping its neighbors. Doing this once produces +# skewed results, so the algorithm is run 4 times (top-down/bottom-up × +# left-biased/right-biased, via reversing the graph's edges and/or each +# rank's order) and the 4 resulting x-coordinates per vertex are combined by +# averaging the two middle (median) values, which Brandes & Köpf show is +# both order- and separation-preserving. +# +function reset_alignment!(g::SugiGraph, layers::Vector{Vector{Int}}) + for (r, layer) in enumerate(layers) + for (p, v) in enumerate(layer) + vert = g.verts[v] + vert.rank = r + vert.pos = p + vert.shift = Inf + vert.align = v + vert.root = v + vert.sink = v + end + end +end + +function is_incident_to_inner_segment(g::SugiGraph, id::Int) + g.verts[id].is_dummy || return false + return any(g.verts[t].is_dummy for t in in_neighbors(g, id)) +end + +function get_inner_segment_upper_neighbor(g::SugiGraph, id::Int) + is_incident_to_inner_segment(g, id) || return nothing + for t in in_neighbors(g, id) + return t + end + return nothing +end + +# Must run after `reset_alignment!` has populated `.pos` for `layers` +# (rust-port #26). The inner loop includes `l_1` itself, not just +# everything strictly before it, or the last vertex of each boundary +# segment never gets checked (rust-port #27). +function mark_type1_conflicts!(g::SugiGraph, layers::Vector{Vector{Int}}) + for r in 1:(length(layers) - 1) + level = layers[r] + next_level = layers[r + 1] + left_dummy_index = 0 + l = 0 + for l_1 in eachindex(next_level) + dummy_candidate = next_level[l_1] + upn = get_inner_segment_upper_neighbor(g, dummy_candidate) + if upn !== nothing + right_dummy_index = g.verts[upn].pos + elseif l_1 == length(next_level) + right_dummy_index = length(level) + else + continue + end + while l < l_1 + vertex = next_level[l + 1] + upper_neighbors = sort!(collect(in_neighbors(g, vertex)); by=u -> g.verts[u].pos) + for un in upper_neighbors + vertex_index = g.verts[un].pos + if vertex_index < left_dummy_index || vertex_index > right_dummy_index + eid = find_edge(g, un, vertex) + g.edges[eid].has_type1_conflict = true + end + end + l += 1 + end + left_dummy_index = right_dummy_index + end + end +end + +function create_vertical_alignments!(g::SugiGraph, layers::Vector{Vector{Int}}) + for layer in layers + r = 0 # sentinel: no neighbor aligned yet + for v in layer + edges = [(eid, g.edges[eid].tail) for eid in g.inn[v] if slack(g, eid, 1) == 0] + isempty(edges) && continue + sort!(edges; by=x -> g.verts[x[2]].pos) + + d = (length(edges) + 1) / 2 - 1 + for m in (floor(Int, d), ceil(Int, d)) + g.verts[v].align == v || continue + eid, median_neighbor = edges[m + 1] + if !g.edges[eid].has_type1_conflict && r < g.verts[median_neighbor].pos + g.verts[median_neighbor].align = v + g.verts[v].root = g.verts[median_neighbor].root + g.verts[v].align = g.verts[v].root + r = g.verts[median_neighbor].pos + end + end + end + end +end + +function compute_block_max_vertex_widths!(g::SugiGraph) + for root in 1:_nv(g) + g.verts[root].root == root || continue + maxw = g.verts[root].width + cur = g.verts[root].align + while cur != root + maxw = max(maxw, g.verts[cur].width) + cur = g.verts[cur].align + end + g.verts[root].block_width = maxw + cur = g.verts[root].align + while cur != root + g.verts[cur].block_width = maxw + cur = g.verts[cur].align + end + end +end + +"""The vertex immediately to the left of `v` within its own rank.""" +pred(g::SugiGraph, v::Int, layers::Vector{Vector{Int}}) = layers[g.verts[v].rank][g.verts[v].pos - 1] + +function place_block!(g::SugiGraph, layers::Vector{Vector{Int}}, root::Int, x::Dict{Int,Float64}) + haskey(x, root) && return + x[root] = 0.0 + w = root + while true + if g.verts[w].pos > 1 + u = g.verts[pred(g, w, layers)].root + place_block!(g, layers, u, x) + g.verts[root].sink == root && (g.verts[root].sink = g.verts[u].sink) + if g.verts[root].sink == g.verts[u].sink + gap = (g.verts[root].block_width + g.verts[u].block_width) * 0.5 + x[root] = max(x[root], x[u] + gap) + end + end + w = g.verts[w].align + w == root && break + end + while g.verts[w].align != root + w = g.verts[w].align + x[w] = x[root] + g.verts[w].sink = g.verts[root].sink + end +end + +function place_blocks(g::SugiGraph, layers::Vector{Vector{Int}}) + x = Dict{Int,Float64}() + for root in 1:_nv(g) + g.verts[root].root == root && place_block!(g, layers, root, x) + end + return x +end + +function do_horizontal_compaction!(g::SugiGraph, layers::Vector{Vector{Int}}) + compute_block_max_vertex_widths!(g) + x = place_blocks(g, layers) + + nlayers = length(layers) + for i in 1:nlayers + v = layers[i][1] + g.verts[v].sink == v || continue + vsink = g.verts[v].sink + g.verts[vsink].shift == Inf && (g.verts[vsink].shift = 0.0) + + j, k = i, 1 + while true + v = layers[j][k] + while g.verts[v].align != g.verts[v].root + v = g.verts[v].align + j += 1 + if g.verts[v].pos > 1 + u = pred(g, v, layers) + gap = (g.verts[v].block_width + g.verts[u].block_width) * 0.5 + distance_v_u = x[v] - (x[u] + gap) + u_sink = g.verts[u].sink + g.verts[u_sink].shift = min(g.verts[u_sink].shift, + g.verts[g.verts[v].sink].shift + distance_v_u) + end + end + k = g.verts[v].pos + 1 + (k > length(layers[j]) || g.verts[v].sink != g.verts[layers[j][k]].sink) && break + end + end + + for v in 1:_nv(g) + x[v] = x[v] + g.verts[g.verts[v].sink].shift + end + return x +end + +""" + create_layouts(g, layers) + +Run the 4-direction Brandes & Köpf alignment and return the resulting list +of 4 coordinate maps (`Dict{Int,Float64}`, vertex id => x coordinate). +""" +function create_layouts(g::SugiGraph, layers::Vector{Vector{Int}}) + reset_alignment!(g, layers) + mark_type1_conflicts!(g, layers) + + layouts = Dict{Int,Float64}[] + cur_layers = [copy(l) for l in layers] + for _vdir in 1:2 + for hdir in 1:2 + reset_alignment!(g, cur_layers) + create_vertical_alignments!(g, cur_layers) + layout = do_horizontal_compaction!(g, cur_layers) + if hdir == 2 # :left + for k in keys(layout) + layout[k] = -layout[k] + end + end + push!(layouts, layout) + for row in cur_layers + reverse!(row) + end + end + reverse_graph!(g) + reverse!(cur_layers) + end + reset_alignment!(g, layers) + return layouts +end + +function align_to_smallest_width_layout!(layouts::Vector{Dict{Int,Float64}}) + isempty(layouts) && return layouts + min_max = map(layouts) do c + lo, hi = extrema(values(c)) + (lo, hi, hi - lo) + end + _, min_width = findmin(x -> x[3], min_max) + + for (i, layout) in enumerate(layouts) + shift = isodd(i) ? min_max[i][1] - min_max[min_width][1] : + min_max[min_width][2] - min_max[i][2] + for k in keys(layout) + layout[k] += shift + end + end + return layouts +end + +"""Average the two median x-coordinates (of the 4 alignment directions) for +each vertex — "the average median is both order and separation preserving" +(Brandes & Köpf, 2002).""" +function calculate_relative_coords(layouts::Vector{Dict{Int,Float64}}) + coords = Dict{Int,Float64}() + for k in keys(layouts[1]) + v = sort!([layouts[1][k], layouts[2][k], layouts[3][k], layouts[4][k]]) + coords[k] = (v[2] + v[3]) / 2.0 + end + return coords +end diff --git a/src/sugiyama/cycles.jl b/src/sugiyama/cycles.jl new file mode 100644 index 0000000..5dfba7a --- /dev/null +++ b/src/sugiyama/cycles.jl @@ -0,0 +1,41 @@ +# Phase 0: cycle removal. +# +# Ranking (phase 1) needs a topological order, so by the time it runs the +# graph must have no directed cycles. We make the graph acyclic by running a +# DFS and reversing every "back edge" found (an edge to a vertex that is +# currently on the DFS stack). + +function remove_cycles!(g::SugiGraph) + n = _nv(g) + visited = falses(n) + onstack = falses(n) + to_reverse = Int[] + + function dfs(v::Int) + visited[v] = true + onstack[v] = true + for eid in copy(g.out[v]) + e = g.edges[eid] + e === nothing && continue + w = e.head + if !visited[w] + dfs(w) + elseif onstack[w] + push!(to_reverse, eid) + end + end + onstack[v] = false + end + + for v in 1:n + visited[v] || dfs(v) + end + + for eid in to_reverse + e = g.edges[eid] + tail, head, weight = e.tail, e.head, e.weight + _rem_edge!(g, eid) + _add_edge!(g, head, tail; weight) + end + return g +end diff --git a/src/sugiyama/graph.jl b/src/sugiyama/graph.jl new file mode 100644 index 0000000..d8bb3e5 --- /dev/null +++ b/src/sugiyama/graph.jl @@ -0,0 +1,131 @@ +# Internal mutable directed graph used to implement the Sugiyama-style +# layered layout (phases 0-3 below). Not a general purpose graph type: it +# only supports what the algorithm needs (add/remove edges and vertices, +# look up edges by endpoints, walk incident edges). Edge ids stay valid +# after removal (removed edges are tombstoned to `nothing` instead of +# compacted), since e.g. network simplex keeps edge ids around across +# mutations. + +mutable struct SugiVertex + rank::Int + pos::Int # position within its layer + low::Int + lim::Int + parent::Union{Nothing,Int} # `nothing` at the root of the feasible tree + is_tree_vertex::Bool + is_dummy::Bool + root::Int + align::Int + shift::Float64 + sink::Int + block_width::Float64 + width::Float64 + height::Float64 +end + +function SugiVertex(; width=1.0, height=1.0, is_dummy=false) + # root/align/sink are placeholders, overwritten with the vertex's own id + # right after construction (see `_add_vertex!`). + SugiVertex(0, 0, 0, 0, nothing, false, is_dummy, 0, 0, Inf, 0, 0.0, width, height) +end + +mutable struct SugiEdge + tail::Int + head::Int + weight::Int + cut_value::Union{Nothing,Int} + is_tree_edge::Bool + has_type1_conflict::Bool +end + +SugiEdge(tail, head; weight=1) = SugiEdge(tail, head, weight, nothing, false, false) + +struct SugiGraph + verts::Vector{SugiVertex} + edges::Vector{Union{Nothing,SugiEdge}} + out::Vector{Vector{Int}} # out[v] = ids of edges with tail == v + inn::Vector{Vector{Int}} # inn[v] = ids of edges with head == v +end + +SugiGraph() = SugiGraph(SugiVertex[], Union{Nothing,SugiEdge}[], Vector{Int}[], Vector{Int}[]) + +_nv(g::SugiGraph) = length(g.verts) + +function _add_vertex!(g::SugiGraph; kwargs...) + push!(g.verts, SugiVertex(; kwargs...)) + push!(g.out, Int[]) + push!(g.inn, Int[]) + id = length(g.verts) + v = g.verts[id] + v.root = id + v.align = id + v.sink = id + return id +end + +function _add_edge!(g::SugiGraph, tail::Int, head::Int; weight=1) + push!(g.edges, SugiEdge(tail, head; weight)) + eid = length(g.edges) + push!(g.out[tail], eid) + push!(g.inn[head], eid) + return eid +end + +function _rem_edge!(g::SugiGraph, eid::Int) + e = g.edges[eid] + g.edges[eid] = nothing + deleteat!(g.out[e.tail], findfirst(==(eid), g.out[e.tail])) + deleteat!(g.inn[e.head], findfirst(==(eid), g.inn[e.head])) + return e +end + +"""Swap the direction of every edge in the graph (in place).""" +function reverse_graph!(g::SugiGraph) + for e in g.edges + e === nothing && continue + e.tail, e.head = e.head, e.tail + end + for v in 1:_nv(g) + g.out[v], g.inn[v] = g.inn[v], g.out[v] + end + return g +end + +out_neighbors(g::SugiGraph, v::Int) = (g.edges[eid].head for eid in g.out[v]) +in_neighbors(g::SugiGraph, v::Int) = (g.edges[eid].tail for eid in g.inn[v]) + +"""Edge ids and the vertex at the other end, considering both directions.""" +function incident(g::SugiGraph, v::Int) + Iterators.flatten((((eid, g.edges[eid].head) for eid in g.out[v]), + ((eid, g.edges[eid].tail) for eid in g.inn[v]))) +end + +function find_edge(g::SugiGraph, tail::Int, head::Int) + for eid in g.out[tail] + g.edges[eid].head == head && return eid + end + return nothing +end + +function find_edge_undirected(g::SugiGraph, a::Int, b::Int) + e = find_edge(g, a, b) + e !== nothing && return e + return find_edge(g, b, a) +end + +function tree_degree(g::SugiGraph, v::Int) + c = 0 + for eid in g.out[v] + g.edges[eid].is_tree_edge && (c += 1) + end + for eid in g.inn[v] + g.edges[eid].is_tree_edge && (c += 1) + end + return c +end + +"""slack of edge `eid`: how much longer it is than the minimum rank length.""" +function slack(g::SugiGraph, eid::Int, minimum_length::Int) + e = g.edges[eid] + return g.verts[e.head].rank - g.verts[e.tail].rank - minimum_length +end diff --git a/src/sugiyama/ordering.jl b/src/sugiyama/ordering.jl new file mode 100644 index 0000000..4993fe0 --- /dev/null +++ b/src/sugiyama/ordering.jl @@ -0,0 +1,234 @@ +# Phase 2: insert dummy vertices for edges spanning more than one rank, then +# reorder vertices within each rank to reduce edge crossings. +# +# Crossing minimization uses the classic median/barycenter bilayer-sweep +# heuristic (Sugiyama, Tagawa & Toda 1981 / Gansner et al. 1993): repeatedly +# sweep down then up through the ranks, reordering each rank by the +# median/average position of its already-fixed neighbors in the previous +# sweep direction, optionally followed by a greedy pairwise `transpose` pass +# that swaps adjacent vertices whenever doing so reduces the local crossing +# count. We keep the best ordering found and stop once a few sweeps in a row +# fail to improve on it (finding the true minimum is NP-hard). + +function insert_dummy_vertices!(g::SugiGraph, minimum_length::Int, dummy_width::Float64) + for eid in [eid for eid in eachindex(g.edges) if g.edges[eid] !== nothing] + slack(g, eid, minimum_length) <= 0 && continue + e = g.edges[eid] + tail, head = e.tail, e.head + _rem_edge!(g, eid) + t = tail + for r in (g.verts[tail].rank + 1):(g.verts[head].rank - 1) + d = _add_vertex!(g; width=dummy_width, height=0.0, is_dummy=true) + g.verts[d].rank = r + _add_edge!(g, t, d) + t = d + end + _add_edge!(g, t, head) + end +end + +# ---- vertex ordering within ranks -------------------------------------- + +mutable struct SugiOrder + layers::Vector{Vector{Int}} + positions::Vector{Int} # position within layer, indexed by vertex id +end + +function SugiOrder(layers::Vector{Vector{Int}}, n::Int) + positions = zeros(Int, n) + for layer in layers, (p, v) in enumerate(layer) + positions[v] = p + end + return SugiOrder(layers, positions) +end + +Base.copy(o::SugiOrder) = SugiOrder([copy(l) for l in o.layers], copy(o.positions)) + +function init_order(g::SugiGraph) + n = _nv(g) + n == 0 && return SugiOrder(Vector{Int}[], 0) + max_rank = maximum(v.rank for v in g.verts) + layers = [Int[] for _ in 1:max_rank] + visited = falses(n) + function dfs(v::Int) + visited[v] && return + visited[v] = true + push!(layers[g.verts[v].rank], v) + for n in collect(out_neighbors(g, v)) + dfs(n) + end + end + for v in 1:n + dfs(v) + end + return SugiOrder(layers, n) +end + +function exchange!(order::SugiOrder, r::Int, i::Int, j::Int) + layer = order.layers[r] + order.positions[layer[i]], order.positions[layer[j]] = j, i + layer[i], layer[j] = layer[j], layer[i] + return order +end + +function barycenter(g::SugiGraph, v::Int, move_down::Bool, positions::Vector{Int}) + neighbors = collect(move_down ? in_neighbors(g, v) : out_neighbors(g, v)) + isempty(neighbors) && return Float64(positions[v]) + return sum(positions[n] for n in neighbors) / length(neighbors) +end + +function median(g::SugiGraph, v::Int, move_down::Bool, positions::Vector{Int}) + neighbors = move_down ? in_neighbors(g, v) : out_neighbors(g, v) + adj = sort!([positions[n] for n in neighbors if abs(g.verts[v].rank - g.verts[n].rank) == 1]) + p = length(adj) + p == 0 && return Inf + m = p ÷ 2 + if isodd(p) + return Float64(adj[m + 1]) + elseif p == 2 + return (adj[1] + adj[2]) / 2.0 + else + left = adj[m] - adj[1] + right = adj[p] - adj[m + 1] + return (adj[m] * right + adj[m + 1] * left) / (left + right) + end +end + +function order_layer(g::SugiGraph, move_down::Bool, cur::SugiOrder, cm::F) where {F} + nlayers = length(cur.layers) + new_layers = [copy(l) for l in cur.layers] + positions = copy(cur.positions) + ranks = move_down ? (2:nlayers) : (nlayers - 1):-1:1 + + for r in ranks + layer = new_layers[r] + scores = Dict(v => cm(g, v, move_down, positions) for v in layer) + sort!(layer; by=v -> scores[v]) + for (p, v) in enumerate(layer) + positions[v] = p + end + end + return SugiOrder(new_layers, positions) +end + +function bilayer_cross_count(g::SugiGraph, order::SugiOrder, rank::Int) + north = order.layers[rank] + south = order.layers[rank + 1] + endpoints = Int[] + for v in north + for n in out_neighbors(g, v) + abs(g.verts[v].rank - g.verts[n].rank) == 1 || continue + push!(endpoints, order.positions[n]) + end + end + return count_crossings(endpoints, length(south)) +end + +"""Count inversions of `endpoints` (positions in `1:south_len`) using a +Fenwick tree; equivalent to (but simpler than) the accumulator-tree method +in the reference implementation.""" +function count_crossings(endpoints::Vector{Int}, south_len::Int) + south_len == 0 && return 0 + bit = zeros(Int, south_len) + function bit_update!(i::Int) + while i <= south_len + bit[i] += 1 + i += i & (-i) + end + end + function bit_query(i::Int) + s = 0 + while i > 0 + s += bit[i] + i -= i & (-i) + end + return s + end + cross = 0 + inserted = 0 + for pos in endpoints + cross += inserted - bit_query(pos) + bit_update!(pos) + inserted += 1 + end + return cross +end + +function total_crossings(g::SugiGraph, order::SugiOrder) + isempty(order.layers) && return 0 + return sum(bilayer_cross_count(g, order, r) for r in 1:(length(order.layers) - 1); init=0) +end + +function cross_count_two_vertices(g::SugiGraph, order::SugiOrder, v::Int, w::Int) + crossings = 0 + for dir in (:in, :out) + v_adj = [order.positions[n] for n in (dir === :in ? in_neighbors(g, v) : out_neighbors(g, v))] + w_adj = [order.positions[n] for n in (dir === :in ? in_neighbors(g, w) : out_neighbors(g, w))] + for i in v_adj, j in w_adj + i > j && (crossings += 1) + end + end + return crossings +end + +function transpose!(g::SugiGraph, order::SugiOrder, move_down::Bool) + nlayers = length(order.layers) + improved = true + ranks = move_down ? (1:nlayers) : (nlayers:-1:1) + while improved + improved = false + for r in ranks + layer = order.layers[r] + for i in 1:(length(layer) - 1) + v, w = layer[i], layer[i + 1] + vw = cross_count_two_vertices(g, order, v, w) + wv = cross_count_two_vertices(g, order, w, v) + if vw > wv + improved = true + exchange!(order, r, i, i + 1) + end + end + end + end +end + +function reduce_crossings_bilayer_sweep(g::SugiGraph, order::SugiOrder, cm::F, + do_transpose::Bool) where {F} + length(order.layers) <= 1 && return order + best_crossings = total_crossings(g, order) + best = copy(order) + cur = order + last_best = 0 + i = 0 + while true + cur = order_layer(g, iseven(i), cur, cm) + do_transpose && transpose!(g, cur, iseven(i)) + crossings = total_crossings(g, cur) + if crossings < best_crossings + best_crossings = crossings + best = copy(cur) + last_best = 0 + else + last_best += 1 + end + last_best == 4 && return best + i += 1 + end +end + +""" + ordering(g, crossing_minimization, transpose) + +Return a `Vector{Vector{Int}}` giving, for every rank, the vertices in that +rank ordered left-to-right so as to (heuristically) minimize edge crossings. +`crossing_minimization` is `:barycenter` or `:median`. +""" +function ordering(g::SugiGraph, crossing_minimization::Symbol, do_transpose::Bool) + order = init_order(g) + cm = crossing_minimization === :barycenter ? barycenter : + crossing_minimization === :median ? median : + throw(ArgumentError("Unknown crossing_minimization $(repr(crossing_minimization)), " * + "must be :barycenter or :median")) + order = reduce_crossings_bilayer_sweep(g, order, cm, do_transpose) + return order.layers +end diff --git a/src/sugiyama/ranking.jl b/src/sugiyama/ranking.jl new file mode 100644 index 0000000..2d0da22 --- /dev/null +++ b/src/sugiyama/ranking.jl @@ -0,0 +1,387 @@ +# Phase 1: assign each vertex a rank (layer). +# +# Follows Gansner, Koutsofios, North & Vo (1993), "A Technique for Drawing +# Directed Graphs": build a feasible spanning tree that is *tight* (every +# tree edge has zero slack), then repeatedly find a tree edge with negative +# cut value and swap it for a non-tree edge until all cut values are +# non-negative, which is optimal (minimizes the total weighted edge length). +# +# `low`/`lim` are a postorder DFS numbering of the feasible tree that lets +# `enter_edge` decide in O(1) whether a candidate edge goes from the tail +# component to the head component after removing a tree edge, without +# re-walking the tree. + +function toposort(g::SugiGraph) + n = _nv(g) + indeg = zeros(Int, n) + for e in g.edges + e === nothing && continue + indeg[e.head] += 1 + end + queue = [v for v in 1:n if indeg[v] == 0] + order = Int[] + sizehint!(order, n) + while !isempty(queue) + v = popfirst!(queue) + push!(order, v) + for eid in g.out[v] + e = g.edges[eid] + e === nothing && continue + indeg[e.head] -= 1 + indeg[e.head] == 0 && push!(queue, e.head) + end + end + length(order) == n || throw(ArgumentError("graph must be acyclic to rank it")) + return order +end + +function init_rank!(g::SugiGraph, minimum_length::Int) + for v in toposort(g) + maxr = nothing + for eid in g.inn[v] + cand = g.verts[g.edges[eid].tail].rank + minimum_length + maxr = maxr === nothing ? cand : max(maxr, cand) + end + maxr !== nothing && (g.verts[v].rank = maxr) + end +end + +function move_vertices_up!(g::SugiGraph, minimum_length::Int) + for v in 1:_nv(g) + maxr = nothing + for eid in g.inn[v] + cand = g.verts[g.edges[eid].tail].rank + minimum_length + maxr = maxr === nothing ? cand : max(maxr, cand) + end + g.verts[v].rank = maxr === nothing ? 0 : maxr + end +end + +function move_vertices_down!(g::SugiGraph, minimum_length::Int) + _nv(g) == 0 && return + max_rank = maximum(v.rank for v in g.verts) + for v in 1:_nv(g) + minr = nothing + for eid in g.out[v] + cand = g.verts[g.edges[eid].head].rank - minimum_length + minr = minr === nothing ? cand : min(minr, cand) + end + g.verts[v].rank = minr === nothing ? max_rank : minr + end +end + +# ---- feasible tree construction ------------------------------------------- + +function tight_tree!(g::SugiGraph, v::Int, visited_edges::Set{Int}, minimum_length::Int) + node_count = 1 + g.verts[v].is_tree_vertex = true + for (eid, other) in collect(incident(g, v)) + eid in visited_edges && continue + push!(visited_edges, eid) + e = g.edges[eid] + if e.is_tree_edge + node_count += tight_tree!(g, other, visited_edges, minimum_length) + elseif slack(g, eid, minimum_length) == 0 && !g.verts[other].is_tree_vertex + e.is_tree_edge = true + node_count += tight_tree!(g, other, visited_edges, minimum_length) + end + end + return node_count +end + +function is_incident_edge(g::SugiGraph, eid::Int) + e = g.edges[eid] + g.verts[e.tail].is_tree_vertex ⊻ g.verts[e.head].is_tree_vertex +end + +function find_non_tight_edge(g::SugiGraph, minimum_length::Int) + best, bestslack = nothing, nothing + for eid in eachindex(g.edges) + e = g.edges[eid] + (e === nothing || e.is_tree_edge) && continue + is_incident_edge(g, eid) || continue + s = slack(g, eid, minimum_length) + if bestslack === nothing || s < bestslack + bestslack, best = s, eid + end + end + best === nothing && throw(ArgumentError("no non-tight edge found while building feasible tree")) + return best +end + +function tighten_edge!(g::SugiGraph, delta::Int) + for v in g.verts + v.is_tree_vertex && (v.rank += delta) + end +end + +function feasible_tree!(g::SugiGraph, minimum_length::Int) + root = 1 + while tight_tree!(g, root, Set{Int}(), minimum_length) < _nv(g) + eid = find_non_tight_edge(g, minimum_length) + e = g.edges[eid] + delta = slack(g, eid, minimum_length) + g.verts[e.head].is_tree_vertex && (delta = -delta) + tighten_edge!(g, delta) + end + init_cutvalues!(g) + init_low_lim!(g) +end + +# ---- cut values ------------------------------------------------------- + +struct NeighborhoodInfo + cut_value_sum::Int + tree_edge_weight_sum::Int + non_tree_edge_weight_sum::Int + missing_v::Union{Nothing,Int} +end + +function get_neighborhood_info(g::SugiGraph, v::Int, dir::Symbol) + cut_value_sum = 0 + tree_edge_weight_sum = 0 + non_tree_edge_weight_sum = 0 + missing_v = nothing + edge_ids = dir === :in ? g.inn[v] : g.out[v] + for eid in edge_ids + e = g.edges[eid] + if !e.is_tree_edge + non_tree_edge_weight_sum += e.weight + elseif e.cut_value !== nothing + cut_value_sum += e.cut_value + tree_edge_weight_sum += e.weight + elseif missing_v === nothing + missing_v = dir === :in ? e.tail : e.head + else + return nothing + end + end + return NeighborhoodInfo(cut_value_sum, tree_edge_weight_sum, non_tree_edge_weight_sum, missing_v) +end + +function calculate_cut_value(edge_weight::Int, incoming::NeighborhoodInfo, outgoing::NeighborhoodInfo) + edge_weight + incoming.non_tree_edge_weight_sum - incoming.cut_value_sum + + incoming.tree_edge_weight_sum - outgoing.non_tree_edge_weight_sum + + outgoing.cut_value_sum - outgoing.tree_edge_weight_sum +end + +function calculate_cut_values!(g::SugiGraph, queue::Vector{Int}) + while !isempty(queue) + v = popfirst!(queue) + incoming = get_neighborhood_info(g, v, :in) + outgoing = get_neighborhood_info(g, v, :out) + (incoming === nothing || outgoing === nothing) && continue + + missing_v = if incoming.missing_v !== nothing && outgoing.missing_v === nothing + incoming.missing_v + elseif incoming.missing_v === nothing && outgoing.missing_v !== nothing + outgoing.missing_v + else + continue + end + + eid = find_edge(g, v, missing_v) + if eid !== nothing + incoming, outgoing = outgoing, incoming + else + eid = find_edge(g, missing_v, v) + end + g.edges[eid].cut_value = calculate_cut_value(g.edges[eid].weight, incoming, outgoing) + push!(queue, missing_v) + end +end + +function leaves(g::SugiGraph) + return [v for v in 1:_nv(g) if tree_degree(g, v) == 1] +end + +init_cutvalues!(g::SugiGraph) = calculate_cut_values!(g, leaves(g)) + +function update_cutvalues!(g::SugiGraph, removed_edge::Int, swap_edge::Int) + lca = remove_outdated_cut_values!(g, swap_edge, removed_edge) + calculate_cut_values!(g, [g.edges[removed_edge].tail]) + return lca +end + +function remove_outdated_cut_values!(g::SugiGraph, swap_edge::Int, removed_edge::Int) + g.edges[removed_edge].cut_value = nothing + w, x = g.edges[swap_edge].tail, g.edges[swap_edge].head + g.verts[w].lim > g.verts[x].lim && ((w, x) = (x, w)) + + lca = clear_cut_values_to_lca!(g, w, x) + + l = x + while l != lca + parent = g.verts[l].parent + eid = find_edge_undirected(g, l, parent) + g.edges[eid].cut_value = nothing + l = parent + end + return lca +end + +"""Clear cut values along the tree path from `w` towards the root, stopping +at (and returning) the least common ancestor of `w` and `x`.""" +function clear_cut_values_to_lca!(g::SugiGraph, w::Int, x::Int) + parent = g.verts[w].parent + parent === nothing && return w + l = w + while true + eid = find_edge_undirected(g, l, parent) + g.edges[eid].cut_value = nothing + l = parent + if (g.verts[l].low <= g.verts[w].lim && g.verts[x].lim <= g.verts[l].lim) || g.verts[l].parent === nothing + return l + end + parent = g.verts[l].parent + end +end + +# ---- low/lim (postorder DFS numbering of the feasible tree) ----------- + +function dfs_low_lim!(g::SugiGraph, v::Int, parent::Union{Nothing,Int}, max_lim::Base.RefValue{Int}, + visited::BitVector) + visited[v] = true + g.verts[v].lim = max_lim[] + g.verts[v].parent = parent + for (eid, other) in collect(incident(g, v)) + if !visited[other] && g.edges[eid].is_tree_edge + max_lim[] -= 1 + dfs_low_lim!(g, other, v, max_lim, visited) + end + end + g.verts[v].low = max_lim[] +end + +function init_low_lim!(g::SugiGraph) + root = 1 + dfs_low_lim!(g, root, nothing, Ref(_nv(g)), falses(_nv(g))) +end + +function update_low_lim!(g::SugiGraph, lca::Int) + parent = g.verts[lca].parent + visited = falses(_nv(g)) + parent !== nothing && (visited[parent] = true) + dfs_low_lim!(g, lca, parent, Ref(g.verts[lca].lim), visited) +end + +# ---- network simplex iteration ----------------------------------------- + +function leave_edge(g::SugiGraph) + for eid in eachindex(g.edges) + e = g.edges[eid] + e === nothing && continue + e.is_tree_edge && e.cut_value !== nothing && e.cut_value < 0 && return eid + end + return nothing +end + +function is_head_to_tail(g::SugiGraph, eid::Int, u::Int, is_root_in_head::Bool) + e = g.edges[eid] + head_in_tail_component = g.verts[u].low <= g.verts[e.head].lim <= g.verts[u].lim + tail_in_head_component = g.verts[u].low <= g.verts[e.tail].lim <= g.verts[u].lim + return (is_root_in_head == head_in_tail_component) && (is_root_in_head != tail_in_head_component) +end + +function enter_edge(g::SugiGraph, eid::Int, minimum_length::Int) + e = g.edges[eid] + u, v = e.tail, e.head + is_root_in_head = g.verts[u].lim < g.verts[v].lim + is_root_in_head || ((u, v) = (v, u)) + + best, bestslack = nothing, nothing + for cid in eachindex(g.edges) + ce = g.edges[cid] + (ce === nothing || ce.is_tree_edge) && continue + is_head_to_tail(g, cid, u, is_root_in_head) || continue + s = slack(g, cid, minimum_length) + if bestslack === nothing || s < bestslack + bestslack, best = s, cid + end + end + best === nothing && throw(ArgumentError("no entering edge found (this should not happen)")) + return best +end + +function exchange!(g::SugiGraph, removed_edge::Int, swap_edge::Int, minimum_length::Int) + g.edges[removed_edge].is_tree_edge = false + g.edges[swap_edge].is_tree_edge = true + lca = update_cutvalues!(g, removed_edge, swap_edge) + update_low_lim!(g, lca) + update_ranks!(g, minimum_length) +end + +function update_neighbor_ranks!(g::SugiGraph, parent::Int, dir::Symbol, coeff::Int, + queue::Vector{Int}, visited::BitVector, minimum_length::Int) + edge_ids = dir === :out ? g.out[parent] : g.inn[parent] + for eid in edge_ids + e = g.edges[eid] + e.is_tree_edge || continue + other = dir === :out ? e.head : e.tail + visited[other] && continue + g.verts[other].rank = g.verts[parent].rank + minimum_length * coeff + push!(queue, other) + visited[other] = true + end +end + +function update_ranks!(g::SugiGraph, minimum_length::Int) + node = 1 + visited = falses(_nv(g)) + visited[node] = true + g.verts[node].rank = 0 + queue = [node] + while !isempty(queue) + parent = popfirst!(queue) + update_neighbor_ranks!(g, parent, :out, 1, queue, visited, minimum_length) + update_neighbor_ranks!(g, parent, :in, -1, queue, visited, minimum_length) + end +end + +"""Shift ranks so the smallest one is `1` (ranks are later used directly as +`layers` row indices, see `coordinates.jl`).""" +function normalize_ranks!(g::SugiGraph) + _nv(g) == 0 && return + shift = 1 - minimum(v.rank for v in g.verts) + for v in g.verts + v.rank += shift + end +end + +function minimize_edge_length!(g::SugiGraph, minimum_length::Int) + feasible_tree!(g, minimum_length) + while (e = leave_edge(g)) !== nothing + swap_edge = enter_edge(g, e, minimum_length) + exchange!(g, e, swap_edge, minimum_length) + end +end + +""" + rank!(g, minimum_length, ranking_type) + +Assign a `rank` to every vertex of `g`. `ranking_type` is one of: +- `:networksimplex`: optimal (minimizes total weighted edge length), via the + network simplex method of Gansner et al. (1993). +- `:longestpath`: move every vertex as far up as possible, then pull sinks + as far down as possible without violating edge lengths. Cheap but tends + to produce wide, uneven layers. +- `:up`: longest path from the sources only. +- `:down`: longest path from the sinks only. +""" +function rank!(g::SugiGraph, minimum_length::Int, ranking_type::Symbol) + init_rank!(g, minimum_length) + if ranking_type === :networksimplex + minimize_edge_length!(g, minimum_length) + elseif ranking_type === :longestpath + move_vertices_up!(g, minimum_length) + move_vertices_down!(g, minimum_length) + elseif ranking_type === :up + move_vertices_up!(g, minimum_length) + elseif ranking_type === :down + move_vertices_down!(g, minimum_length) + else + throw(ArgumentError("Unknown ranking_type $(repr(ranking_type)), must be one of " * + ":networksimplex, :longestpath, :up, :down")) + end + normalize_ranks!(g) +end diff --git a/test/runtests.jl b/test/runtests.jl index bfe74c0..213b21a 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -341,6 +341,8 @@ jagmesh_adj = jagmesh() end end + include("sugiyama_test.jl") + @testset "Testing Square Grid Layout" begin println("SquareGrid") @testset "Testing col length" begin diff --git a/test/sugiyama_test.jl b/test/sugiyama_test.jl new file mode 100644 index 0000000..516124f --- /dev/null +++ b/test/sugiyama_test.jl @@ -0,0 +1,457 @@ +@testset "Testing Sugiyama Layout" begin + using NetworkLayout: SugiGraph, _add_vertex!, _add_edge!, find_edge, _nv, slack + using NetworkLayout: remove_cycles!, rank!, init_cutvalues!, init_low_lim! + using NetworkLayout: leave_edge, enter_edge + using NetworkLayout: insert_dummy_vertices!, bilayer_cross_count, total_crossings, SugiOrder + using NetworkLayout: order_layer, barycenter, ordering + using NetworkLayout: reset_alignment!, mark_type1_conflicts!, create_vertical_alignments!, + reverse_graph!, place_blocks + + println("Sugiyama") + + # ---- helpers: build a SugiGraph from 0-based (tail, head) edges -------- + function sugi_from_edges(edges::Vector{Tuple{Int,Int}}, n::Int) + g = SugiGraph() + for _ in 1:n + _add_vertex!(g) + end + for (t, h) in edges + _add_edge!(g, t + 1, h + 1) + end + return g + end + + function set_tree_edges!(g, tree_edges::Vector{Tuple{Int,Int}}) + for (t, h) in tree_edges + g.verts[t + 1].is_tree_vertex = true + g.verts[h + 1].is_tree_vertex = true + g.edges[find_edge(g, t + 1, h + 1)].is_tree_edge = true + end + end + + @testset "Sugiyama construction" begin + algo = Sugiyama() + @test algo isa Sugiyama{Float64} + @test algo.ranking_type == :networksimplex + @test algo.crossing_minimization == :barycenter + @test algo.direction == :down + + algo = Sugiyama(; Ptype=Float32) + @test algo isa Sugiyama{Float32} + + @test_throws ArgumentError Sugiyama(; direction=:sideways) + @test_throws ArgumentError Sugiyama(; ranking_type=:foo) + @test_throws ArgumentError Sugiyama(; crossing_minimization=:foo) + @test_throws ArgumentError Sugiyama(; minimum_length=0) + end + + # ---- Phase 1: ranking / network simplex, ported from rust-sugiyama's -- + # ---- p1_layering test fixtures (src/algorithm/p1_layering/tests.rs) -- + @testset "ranking: network simplex (ported rust fixtures)" begin + EXAMPLE_GRAPH = [(0, 1), (1, 2), (2, 3), (3, 7), (4, 6), (5, 6), (6, 7), (0, 4), (0, 5)] + FEASIBLE_TREE_NEG = [(0, 1), (1, 2), (2, 3), (3, 7), (4, 6), (5, 6), (6, 7)] + FEASIBLE_TREE_POS = [(0, 1), (0, 4), (1, 2), (2, 3), (3, 7), (4, 6), (5, 6)] + LOW_LIM_GRAPH = [(0, 1), (1, 2), (1, 3), (0, 4), (4, 5), (5, 6), (4, 7), (4, 8)] + + @testset "cut values, tree with a negative cut value" begin + g = sugi_from_edges(EXAMPLE_GRAPH, 8) + set_tree_edges!(g, FEASIBLE_TREE_NEG) + init_cutvalues!(g) + expected = [3, 3, 3, 3, 0, 0, -1] + for (i, (t, h)) in enumerate(FEASIBLE_TREE_NEG) + @test g.edges[find_edge(g, t + 1, h + 1)].cut_value == expected[i] + end + end + + @testset "cut values, all-positive tree" begin + g = sugi_from_edges(EXAMPLE_GRAPH, 8) + set_tree_edges!(g, FEASIBLE_TREE_POS) + init_cutvalues!(g) + expected = [2, 1, 2, 2, 2, 1, 0] + for (i, (t, h)) in enumerate(FEASIBLE_TREE_POS) + @test g.edges[find_edge(g, t + 1, h + 1)].cut_value == expected[i] + end + end + + @testset "low/lim: parent structure and subtree-containment invariant" begin + # the exact low/lim *numbers* are traversal-order dependent (an + # implementation detail); what matters is that they correctly + # encode subtree containment, which is what network simplex relies on. + g = sugi_from_edges(LOW_LIM_GRAPH, 9) + set_tree_edges!(g, LOW_LIM_GRAPH) + init_low_lim!(g) + + expected_parent = Dict(0 => nothing, 1 => 0, 2 => 1, 3 => 1, 4 => 0, 5 => 4, + 6 => 5, 7 => 4, 8 => 4) + for (id, parent) in expected_parent + @test g.verts[id + 1].parent == (parent === nothing ? nothing : parent + 1) + end + + function subtree_of(root, tree_edges, n) + children = [Int[] for _ in 1:n] + for (t, h) in tree_edges + push!(children[t + 1], h + 1) + end + seen, stack = Set([root]), [root] + while !isempty(stack) + v = pop!(stack) + for c in children[v] + c in seen || (push!(seen, c); push!(stack, c)) + end + end + return seen + end + + for x in 0:8, y in 0:8 + sub = subtree_of(x + 1, LOW_LIM_GRAPH, 9) + vx, vy = g.verts[x + 1], g.verts[y + 1] + @test ((y + 1) in sub) == (vx.low <= vy.lim <= vx.lim) + end + end + + @testset "leave_edge / enter_edge" begin + g = sugi_from_edges(EXAMPLE_GRAPH, 8) + set_tree_edges!(g, FEASIBLE_TREE_NEG) + init_cutvalues!(g) + init_low_lim!(g) + + e = leave_edge(g) + @test e !== nothing + @test (g.edges[e].tail, g.edges[e].head) == (7, 8) # rust node 6->7 + + swap = enter_edge(g, e, 1) + @test g.edges[swap].tail == 1 # rust node 0 + @test g.edges[swap].head in (5, 6) # rust node 4 or 5 + + g2 = sugi_from_edges(EXAMPLE_GRAPH, 8) + set_tree_edges!(g2, FEASIBLE_TREE_POS) + init_cutvalues!(g2) + init_low_lim!(g2) + @test leave_edge(g2) === nothing + end + + @testset "full network simplex is optimal & consistent" begin + function is_correct(g, minimum_length) + tree_cv = [g.edges[eid].cut_value for eid in eachindex(g.edges) + if g.edges[eid] !== nothing && g.edges[eid].is_tree_edge] + all(cv -> cv !== nothing && cv >= 0, tree_cv) || return false + tree_slacks = [slack(g, eid, minimum_length) for eid in eachindex(g.edges) + if g.edges[eid] !== nothing && g.edges[eid].is_tree_edge] + all(==(0), tree_slacks) || return false + return minimum(v.rank for v in g.verts) == 1 + end + + g = sugi_from_edges(EXAMPLE_GRAPH, 8) + rank!(g, 1, :networksimplex) + @test is_correct(g, 1) + + # a bigger, denser DAG + import Random + rng = Random.MersenneTwister(42) + n = 80 + edges = Tuple{Int,Int}[] + for i in 0:(n - 2), _ in 1:2 + push!(edges, (i, rand(rng, (i + 1):min(i + 5, n - 1)))) + end + g2 = sugi_from_edges(unique(edges), n) + rank!(g2, 1, :networksimplex) + @test is_correct(g2, 1) + end + + @testset "verify_looks_good graph matches rust's expected width/height" begin + # rust's `verify_looks_good` test asserts width==4.0, height==6.0 + # (max nodes on a rank, and number of ranks) for this graph. + edges = [(0, 1), (1, 2), (2, 3), (2, 4), (3, 5), (3, 6), (3, 7), (3, 8), + (4, 5), (4, 6), (4, 7), (4, 8), (5, 9), (6, 9), (7, 9), (8, 9)] + g = sugi_from_edges(edges, 10) + remove_cycles!(g) + rank!(g, 1, :networksimplex) + @test length(unique(v.rank for v in g.verts)) == 6 + insert_dummy_vertices!(g, 1, 0.0) + layers = ordering(g, :barycenter, true) + @test maximum(length.(layers)) == 4 + end + end + + # ---- Phase 2: crossing reduction, ported from p2_reduce_crossings/tests.rs + @testset "ordering / crossing minimization (ported rust fixtures)" begin + @testset "bilayer_cross_count" begin + g = sugi_from_edges([(0, 4), (1, 3), (2, 3)], 5) # 0-based ids -> n0,n1,n2 / s0,s1 + for v in 1:3 + g.verts[v].rank = 0 + end + for v in 4:5 + g.verts[v].rank = 1 + end + order = SugiOrder([[1, 2, 3], [4, 5]], 5) + @test bilayer_cross_count(g, order, 1) == 2 + + g2 = sugi_from_edges([(0, 7), (1, 6), (2, 5), (3, 4)], 8) + for v in 1:4 + g2.verts[v].rank = 0 + end + for v in 5:8 + g2.verts[v].rank = 1 + end + order2 = SugiOrder([[1, 2, 3, 4], [5, 6, 7, 8]], 8) + @test bilayer_cross_count(g2, order2, 1) == 6 + end + + @testset "total crossings" begin + edges = [(0, 6), (1, 7), (1, 8), (2, 6), (2, 9), (2, 10), (3, 6), (3, 9), + (4, 9), (5, 8), (5, 10)] + g = sugi_from_edges(edges, 11) + for v in 1:6 + g.verts[v].rank = 0 + end + for v in 7:11 + g.verts[v].rank = 1 + end + order = SugiOrder([collect(1:6), collect(7:11)], 11) + @test total_crossings(g, order) == 12 + end + + @testset "insert_dummy_vertices" begin + edges = [(0, 1), (1, 2), (2, 3), (3, 7), (4, 6), (5, 6), (6, 7), (0, 4), (0, 5)] + g = sugi_from_edges(edges, 8) + for (v, r) in [(0, 0), (1, 1), (2, 2), (3, 3), (4, 1), (5, 1), (6, 2), (7, 4)] + g.verts[v + 1].rank = r + end + insert_dummy_vertices!(g, 1, 0.0) + @test _nv(g) == 9 + @test count(v -> v.is_dummy, g.verts) == 1 + end + + @testset "barycenter reordering" begin + g = sugi_from_edges([(0, 8), (1, 8), (2, 8), (2, 10), (3, 9), (4, 10), (5, 10), + (6, 11), (7, 12)], 13) + for v in 1:8 + g.verts[v].rank = 0 + end + for v in 9:13 + g.verts[v].rank = 1 + end + inner = [[1, 3, 5, 4, 7, 8, 2, 6], collect(9:13)] + order = SugiOrder(inner, 13) + result = order_layer(g, false, order, barycenter) + @test result.layers[1] == collect(1:8) + end + end + + # ---- Phase 3: Brandes & Köpf alignment, ported from + # ---- p3_calculate_coordinates/tests.rs's `create_test_layout` fixture + @testset "coordinate assignment (ported rust fixtures)" begin + edges30 = [(0, 2), (0, 6), (0, 18), (1, 16), (1, 17), (3, 8), (16, 8), (4, 8), + (17, 19), (18, 20), (5, 8), (5, 9), (6, 8), (6, 21), (7, 10), (7, 11), + (7, 12), (19, 23), (20, 24), (21, 12), (9, 22), (9, 25), (10, 13), + (10, 14), (11, 14), (22, 13), (23, 15), (24, 15), (12, 15), (25, 15)] + layers0 = [[0, 1], [2, 3, 16, 4, 17, 18, 5, 6], [7, 8, 19, 20, 21, 9], + [10, 11, 22, 23, 24, 12, 25], [13, 14, 15]] + + function fixture() + g = sugi_from_edges(edges30, 26) + layers = [[v + 1 for v in row] for row in layers0] + for (rank, row) in enumerate(layers), (pos, v) in enumerate(row) + vert = g.verts[v] + vert.rank, vert.pos = rank, pos + vert.root = vert.align = vert.sink = v + vert.is_dummy = (v - 1) >= 16 + vert.width = vert.height = (v - 1) >= 16 ? 1.0 : 10.0 + end + return g, layers + end + root_of(g, id) = g.verts[id + 1].root - 1 + + @testset "type-1 conflicts" begin + g, l = fixture() + mark_type1_conflicts!(g, l) + for (t, h) in [(6, 8), (7, 12), (5, 8), (9, 22)] + @test g.edges[find_edge(g, t + 1, h + 1)].has_type1_conflict + end + end + + @testset "type-1 conflicts: last vertex of a row must be checked (rust-port #27)" begin + # Two ranks: upper = [A,B,C,D,Dm] (Dm a dummy continuing a chain into + # the lower rank), lower = [E,Dm2,F] (Dm2 continues that chain, F is + # an ordinary vertex and also the last position in the lower rank). + # A->F should be flagged: drawing it straight would cross the Dm-Dm2 + # inner segment. + g = SugiGraph() + for _ in 1:8 + _add_vertex!(g) + end + A, B, C, D, Dm, E, Dm2, F = 1:8 + _add_edge!(g, B, E) # ordinary edge, well inside any window + _add_edge!(g, Dm, Dm2) # inner segment (dummy -> dummy) + _add_edge!(g, A, F) # crosses the inner segment: must be flagged + + layers = [[A, B, C, D, Dm], [E, Dm2, F]] + g.verts[Dm].is_dummy = true + g.verts[Dm2].is_dummy = true + reset_alignment!(g, layers) + mark_type1_conflicts!(g, layers) + + @test !g.edges[find_edge(g, B, E)].has_type1_conflict + @test g.edges[find_edge(g, A, F)].has_type1_conflict + end + + @testset "down-right alignment (exact root & align)" begin + g, l = fixture() + mark_type1_conflicts!(g, l) + reset_alignment!(g, l) + create_vertical_alignments!(g, l) + + exp_root = Dict(0=>0,1=>1,2=>0,3=>3,4=>4,5=>5,6=>6,7=>7,8=>4,9=>9,10=>7, + 11=>11,12=>6,13=>7,14=>11,15=>18,16=>1,17=>17,18=>18,19=>17, + 20=>18,21=>6,22=>22,23=>17,24=>18,25=>9) + for (id, r) in exp_root + @test root_of(g, id) == r + end + + exp_align = Dict(0=>2,1=>16,2=>0,3=>3,4=>8,5=>5,6=>21,7=>10,8=>4,9=>25,10=>13, + 11=>14,12=>6,13=>7,14=>11,15=>18,16=>1,17=>19,18=>20,19=>23, + 20=>24,21=>12,22=>22,23=>17,24=>15,25=>9) + for (id, a) in exp_align + @test g.verts[id + 1].align - 1 == a + end + end + + @testset "down-left / up-right / up-left alignment (block membership)" begin + g, l = fixture() + mark_type1_conflicts!(g, l) + foreach(reverse!, l) + reset_alignment!(g, l) + create_vertical_alignments!(g, l) + for (root, members) in Dict(0=>[0,6], 4=>[4,8], 17=>[17,19,23], 18=>[18,20,24], + 5=>[5,9,25], 7=>[7,11,14], 21=>[21,12,15], 10=>[10,13]) + for m in members + @test root_of(g, m) == root + end + end + + g2, l2 = fixture() + mark_type1_conflicts!(g2, l2) + reverse_graph!(g2) + reverse!(l2) + reset_alignment!(g2, l2) + create_vertical_alignments!(g2, l2) + for (root, members) in Dict(13=>[13,10], 14=>[14,11,7], 15=>[15,23,19,17], + 24=>[24,20,18,0], 12=>[12,21], 25=>[25,9,5], 8=>[8,3]) + for m in members + @test root_of(g2, m) == root + end + end + + g3, l3 = fixture() + mark_type1_conflicts!(g3, l3) + reverse_graph!(g3) + reverse!(l3) + foreach(reverse!, l3) + reset_alignment!(g3, l3) + create_vertical_alignments!(g3, l3) + for (root, members) in Dict(15=>[15,25,9], 13=>[13,22], 12=>[12,21,6], + 24=>[24,20,18], 23=>[23,19,17,1], 11=>[11,7], 8=>[8,4]) + for m in members + @test root_of(g3, m) == root + end + end + end + + @testset "place_blocks sinks" begin + g, l = fixture() + mark_type1_conflicts!(g, l) + create_vertical_alignments!(g, l) + x = place_blocks(g, l) + @test length(x) == 26 + for v in [0,1,2,3,4,5,6,8,9,12,15,16,17,18,19,20,21,23,24,25] + @test g.verts[v + 1].sink - 1 == 0 + end + for v in [7,10,11,22,13,14] + @test g.verts[v + 1].sink - 1 == 7 + end + end + end + + # ---- End-to-end API-level tests ----------------------------------- + @testset "end-to-end layout()" begin + @testset "empty and trivial graphs" begin + @test Sugiyama()(zeros(Int, 0, 0)) == Point{2,Float64}[] + @test length(Sugiyama()(zeros(Int, 1, 1))) == 1 + end + + @testset "self loops are ignored" begin + adj = [1 1; 0 0] + pos = Sugiyama()(adj) + @test length(pos) == 2 + @test all(isfinite, pos[1]) && all(isfinite, pos[2]) + end + + @testset "cyclic graphs are handled (implicit edge reversal)" begin + adj = zeros(Int, 5, 5) + for (i, j) in [(1, 2), (2, 3), (3, 1), (3, 4), (4, 5), (5, 3)] + adj[i, j] = 1 + end + pos = Sugiyama()(adj) + @test length(pos) == 5 + @test all(p -> all(isfinite, p), pos) + end + + @testset "disconnected components are tiled without overlap" begin + adj = [0 1 0 0 0; + 0 0 0 0 0; + 0 0 0 1 0; + 0 0 0 0 1; + 0 0 0 0 0] + pos = Sugiyama()(adj) + @test length(pos) == 5 + @test length(unique(pos)) == 5 + end + + @testset "wheel_graph via Graphs.jl" begin + g = wheel_graph(10) + dirg = SimpleDiGraph(collect(edges(g))) + adj = adjacency_matrix(dirg) + pos = @time Sugiyama()(adj) + @test typeof(pos) == Vector{Point{2,Float64}} + @test pos == sugiyama(adj) + @test pos == Sugiyama()(dirg) + end + + @testset "Ptype and direction keywords" begin + adj = adjacency_matrix(SimpleDiGraph(path_digraph(5))) + pos = Sugiyama(; Ptype=Float32)(adj) + @test typeof(pos) == Vector{Point{2,Float32}} + + for d in (:down, :up, :left, :right) + pos = Sugiyama(; direction=d)(adj) + @test length(pos) == 5 + @test all(p -> all(isfinite, p), pos) + end + # :down and :up should be vertical mirror images (same |y| magnitudes) + pd = Sugiyama(; direction=:down)(adj) + pu = Sugiyama(; direction=:up)(adj) + @test getindex.(pd, 2) == -getindex.(pu, 2) + end + + @testset "ranking_type / crossing_minimization / transpose combinations" begin + adj = adjacency_matrix(SimpleDiGraph(path_digraph(6))) + for rt in (:networksimplex, :longestpath, :up, :down), + cm in (:barycenter, :median), tr in (true, false) + + pos = Sugiyama(; ranking_type=rt, crossing_minimization=cm, transpose=tr)(adj) + @test length(pos) == 6 + end + end + + @testset "nodesize / nodespacing / dummysize keywords" begin + adj = adjacency_matrix(SimpleDiGraph(path_digraph(4))) + pos = Sugiyama(; nodesize=[1.0, 2.0, 0.5], nodespacing=2.0, dummysize=0.3)(adj) + @test length(pos) == 4 + end + end + + @testset "assert square" begin + M1 = rand(2, 4) + @test_throws ArgumentError sugiyama(M1) + end +end