Performance of weighed graphs

I work on combinarotial optimization, so I am looking data structures storing weighed graph.
At https://juliagraphs.org/, I found SimpleWeightedGraphs, MetaGraphs, MetaGraphsNext.
However, neither of them seems to me to be reasonably fast. For a small demonstration,
consider the following comparison on Bellman-Ford algorithm.

using Graphs
using SimpleWeightedGraphs
using MetaGraphs
using MetaGraphsNext
using BenchmarkTools

struct DataGraph <: AbstractGraph{Int}
    fadjlist::Vector{Vector{Tuple{Int,Int}}}
end
DataGraph(nvg::Int) = DataGraph([ Vector{Tuple{Int,Int}}() for _ in 1:nvg ])

Graphs.nv(g::DataGraph) = length(g.fadjlist)

Graphs.vertices(g::DataGraph) = 1:nv(g)

incident(g::DataGraph, u::Int) = g.fadjlist[u]

function Graphs.add_edge!(g::DataGraph, u::Int, v::Int, w::Int)
    push!(g.fadjlist[u], (v,w))
    push!(g.fadjlist[v], (u,w))
end

# Source: https://github.com/JuliaGraphs/Graphs.jl/blob/master/src/shortestpaths/bellman-ford.jl
function Graphs.bellman_ford_shortest_paths(graph::DataGraph, sources::AbstractVector{<:Integer})
    nvg = nv(graph)
    active = falses(nvg)
    active[sources] .= true
    dists = fill(typemax(Int), nvg)
    parents = zeros(Int, nvg)
    dists[sources] .= 0
    no_changes = false
    new_active = falses(nvg)

    for i in vertices(graph)
        no_changes = true
        new_active .= false
        for u in vertices(graph)[active]
            for (v,w) in incident(graph, u) # Changed
                relax_dist = w + dists[u]   # Changed
                if dists[v] > relax_dist
                    dists[v] = relax_dist
                    parents[v] = u
                    no_changes = false
                    new_active[v] = true
                end
            end
        end
        if no_changes
            break
        end
        active, new_active = new_active, active
    end
    no_changes || throw(NegativeCycleError())
    return Graphs.BellmanFordState(parents, dists)
end

function main()
    println("Loading data")

#   https://data.4tu.nl/articles/dataset/Network_instance_Chordal_fixed_treewidth_edition_1/12697523/1
    filename = "4tu_chordal/chordal-varNodes-3125,211,1.dimacs"; srcs = [3125] # 3125 & 637009

    all_edges = [tuple(parse.(Int, split(line)[2:4])...) for line in eachline(filename) if line[1] == 'a' ]

    println("Creating graphs")
    swg = SimpleWeightedGraph([ x[1] for x in all_edges ], [ x[2] for x in all_edges ], [ x[3] for x in all_edges ]; combine = max)
    println("Vertices $(nv(swg)), original edges $(length(all_edges)), half $(length(all_edges)/2), unique $(ne(swg))")
    empty!(all_edges)

    dg = DataGraph(nv(swg))
    mg = MetaGraphs.MetaGraph(nv(swg))
    mgn = MetaGraphsNext.MetaGraph(Graph(), Label = Int64, VertexData = Int64, EdgeData = Int64, weight_function = identity)
    for i in 1:nv(swg)
        add_vertex!(mgn, i, i)
    end
    for e in edges(swg)
        u,v = src(e),dst(e)
        w = get_weight(swg, u, v)
        add_edge!(dg, u, v, w)
        add_edge!(mg, u, v, :weight, w)
        add_edge!(mgn, u, v, w)
    end

    println("Precompiling")
    state1 = bellman_ford_shortest_paths(dg, srcs)
    state2 = bellman_ford_shortest_paths(swg, srcs)
    @assert all(state1.dists .== state2.dists)
    state2 = bellman_ford_shortest_paths(mg, srcs)
    @assert all(state1.dists .== state2.dists)
    state2 = bellman_ford_shortest_paths(mgn, srcs)
    @assert all(state1.dists .== state2.dists)
    state1 = state2 = nothing

    println("SimpleWeightedGraph")
    display(@benchmark bellman_ford_shortest_paths($swg, $srcs))
    println()

    println("MetaGraph")
    display(@benchmark bellman_ford_shortest_paths($mg, $srcs))
    println()

    println("MetaGraphNext")
    display(@benchmark bellman_ford_shortest_paths($mgn, $srcs))
    println()

    println("DataGraph")
    display(@benchmark bellman_ford_shortest_paths($dg, $srcs))
    println()

    println("SimpleWeightedGraph")
    display(@benchmark bellman_ford_shortest_paths($swg, $srcs))
    println()
end

main()

The tested graph contains 3125 vertices and 637009 edges. Median running times are

MetaGraph: 1.720 s
MetaGraphNext: 1.097 s
SimpleWeightedGraph: 128.038 ms
My DataGraph: 4.403 ms

The diffence is caused by the way weights are stored: MetaGraph and MetaGraphNext
uses hash tables and SimpleWeightedGraph uses binary search in a sparse matrix.

I know that my graph representation may not be ideal in all situations.
This is just a simple demonstration of the potential for significant improvements.
There are many algorithms accessing all edges incident to a vertex together with its weight or other data.
In order to access the weight of a given edge, this example can be extended e.g. by sorting the list of neighbors or adding a hash table.

Are there some other packages for weighed graphs? Do I miss something?

Tested on
Julia 1.7.0
[86223c79] Graphs v1.7.1
[626554b9] MetaGraphs v0.7.1
[fa8bd995] MetaGraphsNext v0.3.0
[47aef6b3] SimpleWeightedGraphs v1.2.1

Hi @Jirka_Fink and thank you for this benchmark!

First a small remark. In order to avoid extending the functions from Graphs.jl to your own types, all you have to do is make sure they satisfy the AbstractGraph interface.
However, it made sense to re-implement the Bellman-Ford algorithms in this case. The shortest path algorithms from Graphs.jl all use an optional argument distmx = weights(g) which must be an AbstractMatrix, and given the structure of DataGraph it would have been wasteful to construct this matrix in the first place.

Now onto performance. If you want a fair comparison, the only reasonable contender is SimpleWeightedGraphs.jl. Indeed, MetaGraphs.jl and MetaGraphsNext.jl are both meant for a more general use, which means they compromise on speed (as a matter of fact, MetaGraphs.jl is not even type stable). There are other packages for graphs with generic metadata (such as SimpleValueGraphs.jl or MetaDataGraphs.jl), but I doubt they will be able to beat SimpleWeightedGraphs.
I was able to reproduce your results on my computer, and the benchmarking code looks pretty clean to me, so I assume your conclusions are valid. A reassuring note is that SimpleGraph (the unweighted default graph format) was still significantly faster than your custom type, even though of course it’s cheating.

As a maintainer of the Graphs.jl ecosystem, I’d warmly welcome a PR from you on this topic. I think it would make sense to add it to SimpleWeightedGraphs.jl: depending on the application, users may need either adjacency list or adjacency matrix storage, and having both in the same place makes it easier to switch and compare.
Pinging my buddies @etienne_dg and @mbesancon so they can voice their opinion. I’m completely swamped right now so I don’t have a lot of time for Graphs.jl work, but it should get better after the summer.

@gdalle, thank you for you answer.

I already considered how to use AbstractGraph interface but I do not know how.
Notice that I make a simple but crucial change in the code:

    for (v,w) in incident(graph, u)
        relax_dist = w + dists[u]

This is the essential key of my improvement.
Actually, it may be better to change the implementation the Bellman-Ford algorithms instead of introducing new graph structure; see another benchmark.

using Graphs
using SimpleWeightedGraphs
using BenchmarkTools

struct DataGraph <: AbstractGraph{Int}
    fadjlist::Vector{Vector{Tuple{Int,Int}}}
end
DataGraph(nvg::Int) = DataGraph([ Vector{Tuple{Int,Int}}() for _ in 1:nvg ])

Graphs.nv(g::DataGraph) = length(g.fadjlist)

Graphs.vertices(g::DataGraph) = 1:nv(g)

incident(g::DataGraph, u::Int) = g.fadjlist[u]

function Graphs.add_edge!(g::DataGraph, u::Int, v::Int, w::Int)
    push!(g.fadjlist[u], (v,w))
    push!(g.fadjlist[v], (u,w))
end

function incident(g::AbstractSimpleWeightedGraph, v::Integer)
    mat = g.weights
    indices = mat.colptr[v]:mat.colptr[v+1]-1
    return zip(view(mat.rowval, indices), view(mat.nzval, indices))
end

function bellman_ford_shortest_paths_incident(graph::AbstractGraph, sources::AbstractVector{<:Integer})
    nvg = nv(graph)
    active = falses(nvg)
    active[sources] .= true
    dists = fill(typemax(Int), nvg)
    parents = zeros(Int, nvg)
    dists[sources] .= 0
    no_changes = false
    new_active = falses(nvg)

    for i in vertices(graph)
        no_changes = true
        new_active .= false
        for u in vertices(graph)[active]
            for (v,w) in incident(graph, u)
                relax_dist = w + dists[u]
                if dists[v] > relax_dist
                    dists[v] = relax_dist
                    parents[v] = u
                    no_changes = false
                    new_active[v] = true
                end
            end
        end
        if no_changes
            break
        end
        active, new_active = new_active, active
    end
    no_changes || throw(NegativeCycleError())
    return Graphs.BellmanFordState(parents, dists)
end

function main()
    println("Loading data")

#   https://data.4tu.nl/articles/dataset/Network_instance_Chordal_fixed_treewidth_edition_1/12697523/1
    filename = "4tu_chordal/chordal-varNodes-3125,211,1.dimacs"; srcs = [3125] # 3125 vertices and 637009 edges

    @time all_edges = [tuple(parse.(Int, split(line)[2:4])...) for line in eachline(filename) if line[1] == 'a' ]

    println("Creating graphs")
    @time swg = SimpleWeightedGraph([ x[1] for x in all_edges ], [ x[2] for x in all_edges ], [ x[3] for x in all_edges ]; combine = max)
    println("Vertices $(nv(swg)), original edges $(length(all_edges)), half $(length(all_edges)/2), unique $(ne(swg))")
    empty!(all_edges)

    dg = DataGraph(nv(swg))
    @time for e in edges(swg)
        u,v = src(e),dst(e)
        w = get_weight(swg, u, v)
        add_edge!(dg, u, v, w)
    end

    println("Precompiling")
    state1 = bellman_ford_shortest_paths(swg, srcs)
    state2 = bellman_ford_shortest_paths_incident(swg, srcs)
    @assert all(state1.dists .== state2.dists)
    state2 = bellman_ford_shortest_paths_incident(dg, srcs)
    @assert all(state1.dists .== state2.dists)
    state1 = state2 = nothing

    println("SimpleWeightedGraph - original")
    display(@benchmark bellman_ford_shortest_paths($swg, $srcs))
    println()

    println("SimpleWeightedGraph - incident")
    display(@benchmark bellman_ford_shortest_paths_incident($swg, $srcs))
    println()

    println("DataGraph")
    display(@benchmark bellman_ford_shortest_paths_incident($dg, $srcs))
    println()

    println("SimpleWeightedGraph - original")
    display(@benchmark bellman_ford_shortest_paths($swg, $srcs))
    println()
end

main()

SimpleWeightedGraph - original: 128.118 ms
SimpleWeightedGraph - incident: 6.410 ms
DataGraph: 4.261 ms

Before sending a PR, we should decide

  • whether my DataGraph is really so beneficial (more graphs should be evaluated to see the real speed up)
  • whether it would be better to have two implementation of the Bellman-Ford algorithms or modify the current one
  • whether a new interface function incident should be introduced (and it’s name). But weights may not be part of the graph, so the function incident may need both the graph and weights.
1 Like

I always wondered why we are using a sparse matrix to store weights instead of using an adjacency list storage. Glad to see someone questioning it. A drawback will probably a slightly slower access to neighbors, but switching to it would also solve other issues like the fact that zero weights are not supported. Maybe there is a tradeoff decision, and we could have both implementations existing. I will open an issue on SimpleWeightedGraphs.

Edit: An other drawback is that querying weight between two vertices will be slower, and generating the weight matrix will be costly.

On the other side, the idea to iterate on the edges rather than the vertices has been discussed recently (here and here), this could allow this type of optimization.

https://github.com/JuliaGraphs/SimpleWeightedGraphs.jl/issues/25

1 Like

There is also SimpleValueGraphs.jl that uses adjacency lists for storing weights.

In that package I use matrix views of graphs, so that one can avoid the costly generation of weight matrices.

But in general, most of the algorithms in Graphs.jl, that take in a weight matrices can also be speed up, by simply iterating over all outgoing weights of the weightmatrix (i.e. the rows).

In order to compare different representations of weighted graphs, I have more benchmarks. Are you interested to see them? If so, where should I upload them? This place does not seems to me to be right place to put them. It may also be interesting to compare it with API graphs 2.

You can put it on Github.

You can find my benchmarks here: GitHub - jirka-fink/JuliaBenchmarkGraphs: Development and benchmark structures for networks in Julia
Results for Bellman-Ford algorithm are here: Results for Bellman-Ford algorithm · Issue #1 · jirka-fink/JuliaBenchmarkGraphs · GitHub

We already discussed that it is necessary to iterate neighbors together with weights. In this case, sparse matrix used in SimpleWeightedGraph is faster than adjacency lists used in my DataGraph for sparse graphs but it is slower for dense graphs. I am not familiar how iterators are internally implemented.

Please, can somebody explain why SimpleWeightedGraph is faster for sparse graphs but slower for dense graphs? Are there faster ways for used iterations? Recommendation how to make faster structures for weighted graphs are also welcomed.

1 Like

Thanks for the benchmark, it looks very interesting, I will dive into it when I will have some time.