julia系列9：图论包

1. 图论框架对比

这篇文章测试了目前开源的六个图论框架，四个是Python&C++，一个是纯Python（基准），还有一个是Julia的LightGraph.jl，发现Julia的运算速度十分惊人（因为一开始本来没打算测Julia的）。并且与其他框架相比，Julia可以显式控制是否使用并行。

• 这6个框架是：
  graph-tool, v2.31 (Peixoto 2014)
  igraph, v0.8.2 (Csardi and Nepusz 2006)
  networkit, v6.1.0 (Staudt, Sazonovs, and Meyerhenke 2016)
  networkx, v2.4 (Hagberg, Swart, and S Chult 2008)
  SNAP, v5.0.0 (Leskovec and Sosič 2016)
  lightgraphs, v2.0-dev (Seth Bromberger and contributors 2017)

• 测试的5个问题如下：
  loading the data
  single source shortest path
  page rank
  k-core decomposition
  strongly connected components

• 测试数据集为：
  Amazon product co-purchasing network from March 2 2003, 262k nodes, 1.2m edges
  Web graph from Google, 875k nodes, 5.1m edges
  Pokec online social network, 1.6m nodes, 30.6m edges

• 结果：
  lightGraph除了loading比较慢之外，在算法上比其他包要快很多。
  

2. 简单入手

2.1 构建图

1. 逐步添加边

using LightGraphs
using GraphPlot

G₁ = Graph(3) # graph with 3 vertices

# make a triangle
add_edge!(G₁, 1, 2)
add_edge!(G₁, 1, 3)
add_edge!(G₁, 2, 3)

gplot(G₁, nodelabel=1:3)

2. 使用矩阵：

A = [0 1 1
    1 0 1
    1 1 0]

G₂ = Graph(A)

3. 随机生成边：

G = DiGraph(3,3)
gplot(G)

4. 预定义的图，叫generator：

for shape in ["bull","house"]
    G = smallgraph(shape)
    display(gplot(G))
end

G = CliqueGraph(3,4)
gplot(G)

5. 使用edge批量添加

el = Edge.([ (1, 3), (1, 5), (3, 1) ])
SimpleDiGraph(el)

7. 可以用嵌套方式：

# metal plate
G = Graph(4)
add_edge!(G, 1, 2)
add_edge!(G, 1, 3)
add_edge!(G, 2, 4)
add_edge!(G, 3, 4)

# airplane skeleton
skeleton = Graph(11)
add_edge!(skeleton, 1, 2)
add_edge!(skeleton, 2, 3)
add_edge!(skeleton, 3, 4)
add_edge!(skeleton, 4, 5)
add_edge!(skeleton, 3, 6)
add_edge!(skeleton, 3, 7)
add_edge!(skeleton, 3, 8)
add_edge!(skeleton, 3, 9)
add_edge!(skeleton, 9, 10)
add_edge!(skeleton, 9, 11)

gplot(cartesian_product(□, skeleton))

2.2 基础函数

G = smallgraph("house")
# 获取节点
nvertices = nv(G) # number of vertices
nedges = ne(G)    # number of edges
for v in vertices(G)
    println("vertex $v")
end
for e in edges(G)
    u, v = src(e), dst(e)
    println("edge $u - $v")
end
adjacency_matrix(G) # 邻接矩阵
incidence_matrix(G) # 关联矩阵
laplacian_matrix(G)  # 拉普拉斯矩阵
gplot(G, nodelabel=1:nvertices, edgelabel=1:nedges) # 绘图

# 增加/删除点和边
add_vertex!(G)
add_edge!(G, 5, 6)
gplot(G, nodelabel=1:nv(G), edgelabel=1:ne(G))
rem_vertex!(G, 1)
gplot(G, nodelabel=1:nv(G), edgelabel=1:ne(G))

2.3 常用方法

# 是否联通
julia> g = SimpleDiGraph([0 1 0 0 0; 0 0 1 0 0; 1 0 0 1 0; 0 0 0 0 1; 0 0 0 1 0])
{5, 6} directed simple Int64 graph

julia> strongly_connected_components(g)
2-element Array{Array{Int64,1},1}:
 [4, 5]
 [1, 2, 3]

# The attracting components are a subset of the strongly connected components in which the components do not have any leaving edges.
julia> attracting_components(g)
1-element Array{Array{Int64,1},1}:
 [4, 5]

# neighborhood(g, v, d, distmx=weights(g))
# Return a vector of each vertex in g at a geodesic distance less than or equal to d, where distances may be specified by distmx.
julia> g = SimpleDiGraph([0 1 0 0 0; 0 0 1 0 0; 1 0 0 1 0; 0 0 0 0 1; 0 0 0 1 0]);

julia> neighborhood(g, 1, 2)
3-element Array{Int64,1}:
 1
 2
 3
julia> neighborhood_dists(g, 1, 2)
3-element Vector{Tuple{Int64, Int64}}:
 (1, 0)
 (2, 1)
 (3, 2)

# is_cyclic(g)
# Return true if graph g contains a cycle.

# minimum spanning tree
using SimpleWeightedGraphs
g = SimpleWeightedGraph(3)  # or use `SimpleWeightedDiGraph` for directed graphs
add_edge!(g, 1, 2, 0.5)
add_edge!(g, 2, 3, 0.8)
add_edge!(g, 1, 3, 2.0)
kruskal_mst(g)
2-element Vector{SimpleWeightedEdge{Int64, Float64}}:
 Edge 1 => 2 with weight 0.5
 Edge 2 => 3 with weight 0.8

using SimpleWeightedGraphs
g = SimpleWeightedGraph(
        [0 1 0 0 0
        1 0 0 0 0
        0 0 0 1.1 0
        0 0 1.1 0 1
        0 0 0 1 0])
kruskal_mst(g)
3-element Vector{SimpleWeightedEdge{Int64, Float64}}:
 Edge 1 => 2 with weight 1.0
 Edge 4 => 5 with weight 1.0
 Edge 3 => 4 with weight 1.1

dijkstra_shortest_paths(g, [1,5])
Graphs.DijkstraState{Float64, Int64}([0, 1, 4, 5, 0], [0.0, 1.0, 2.1, 1.0, 0.0], [Int64[], Int64[], Int64[], Int64[], Int64[]], [1.0, 1.0, 1.0, 1.0, 1.0], Int64[])

2.4 最大费用流

julia> using Graphs, GraphsFlows
julia> flow_graph = Graphs.DiGraph(8) # Create a flow graph
julia> flow_edges = [
    (1,2,10),(1,3,5),(1,4,15),(2,3,4),(2,5,9),
    (2,6,15),(3,4,4),(3,6,8),(4,7,16),(5,6,15),
    (5,8,10),(6,7,15),(6,8,10),(7,3,6),(7,8,10)
]

julia> capacity_matrix = zeros(Int, 8, 8)  # Create a capacity matrix

julia> for e in flow_edges
    u, v, f = e
    Graphs.add_edge!(flow_graph, u, v)
    capacity_matrix[u,v] = f
end

julia> f, F = maximum_flow(flow_graph, 1, 8) # Run default maximum_flow (push-relabel) without the capacity_matrix

julia> f, F = maximum_flow(flow_graph, 1, 8, capacity_matrix) # Run default maximum_flow with the capacity_matrix

julia> f, F = maximum_flow(flow_graph, 1, 8, capacity_matrix, algorithm=EdmondsKarpAlgorithm()) # Run Edmonds-Karp algorithm

julia> f, F = maximum_flow(flow_graph, 1, 8, capacity_matrix, algorithm=DinicAlgorithm()) # Run Dinic's algorithm

julia> f, F, labels = maximum_flow(flow_graph, 1, 8, capacity_matrix, algorithm=BoykovKolmogorovAlgorithm()) # Run Boykov-Kolmogorov algorithm

julia> mincut(flow_graph, 1, 8, capacity_matrix)