nx.triangles(G)在具有大約 15 萬個節點和 200 萬條邊的無向圖上計算目前非常慢(80 小時的規模)。如果節點度分布高度偏斜,使用下面的程式計算三角形有什么問題嗎?
import networkx as nx
def largest_degree_node(G):
# this was improved using suggestion by Stef in the comments
return max(G.degree(), key=lambda x: x[1])[0]
def count_triangles(G):
G=G.copy()
triangle_counts = 0
while len(G.nodes()):
focal_node = largest_degree_node(G)
triangle_counts = nx.triangles(G, nodes=[focal_node])[focal_node]
G.remove_node(focal_node)
return triangle_counts
G = nx.erdos_renyi_graph(1000, 0.1)
# compute triangles with nx
triangles_nx = int(sum(v for k, v in nx.triangles(G).items()) / 3)
# compute triangles iteratively
triangles_iterative = count_triangles(G)
# assertion passes
assert int(triangles_nx) == int(triangles_iterative)
斷言通過了,但我擔心在某些邊緣情況下這種迭代方法不起作用。
uj5u.com熱心網友回復:
假設圖是無向的(即G.is_directed() == False),可以通過查找既是鄰居的鄰居又是同一節點的直接鄰居的節點來有效地找到三角形的數量。預先計算和預先過濾節點的鄰居,以便每個三角形只計算一次有助于大大提高執行時間。這是代碼:
nodeNeighbours = {
# The filtering of the set ensure each triangle is only computed once
node: set(n for n in edgeInfos.keys() if n > node)
for node, edgeInfos in G.adjacency()
}
triangleCount = sum(
len(neighbours & nodeNeighbours[node2])
for node1, neighbours in nodeNeighbours.items()
for node2 in neighbours
)
上面的代碼比示例圖上的原始迭代解決方案快約12 倍。和最多快72倍的nx.erdos_renyi_graph(15000, 0.005)。
uj5u.com熱心網友回復:
甲更快的解決方案在于使用Numba并與執行計算的多個執行緒。不幸的是,這種方法要復雜得多(因此是這個單獨的答案)。
思路是先把圖轉換成Numpy陣列,然后平衡執行緒之間的作業,最后應用替代答案中提供的演算法(使用多執行緒和Numpy陣列)。這些圖使用四個 Numpy 陣列進行編碼:
nodes: 包含所有節點 ID;allNeighbours: 包含所有節點的所有鄰居;neighbourSlices:包含每個節點的陣列切片(開始 結束索引)allNeighbours(以便能夠獲取節點的鄰居);nodeIdToPos: 包含nodes基于節點ID的節點索引。
這是代碼:
# Split the nodes so the work (based on the slices) is balanced.
# Returns the (start,stop) slices of the nodes
@nb.njit('(int_[:,::1], int_)')
def splitSlices(neighbourSlices, count):
n = np.int64(neighbourSlices.shape[0])
m = neighbourSlices[-1, 1] if neighbourSlices.size > 0 else 0
workSize = np.empty((count, 2), dtype=np.int_)
curPos = 0
for i in range(count):
for j in range(curPos, n):
stop = neighbourSlices[j, 1]
if stop >= m * (i 1) // count:
workSize[i, 0] = curPos
curPos = j 1
workSize[i, 1] = curPos
break
if count > 0:
workSize[0, 0] = 0
workSize[-1, 1] = n
return workSize
@nb.njit('(int_[::1], int_[:,::1], int_[::1])', parallel=True)
def filterNeighbours(nodes, neighbourSlices, allNeighbours):
outNeighbourSlices = np.empty(neighbourSlices.shape, dtype=np.int_)
outAllNeighbours = np.empty(allNeighbours.size//2, dtype=np.int_)
curPos = 0
for i in range(nodes.size):
curNode = nodes[i]
start, stop = neighbourSlices[i]
outNeighbourSlices[i, 0] = curPos
for neighbour in allNeighbours[start:stop]:
if neighbour > curNode:
outAllNeighbours[curPos] = neighbour
curPos = 1
outNeighbourSlices[i, 1] = curPos
threadCount = nb.np.ufunc.parallel.get_num_threads()
nodeSlides = splitSlices(outNeighbourSlices, threadCount)
for threadId in nb.prange(threadCount):
start, stop = nodeSlides[threadId]
for i in range(start, stop):
start, stop = outNeighbourSlices[i]
outAllNeighbours[start:stop].sort()
return (outNeighbourSlices, outAllNeighbours)
@nb.njit('int_(int_[::1], int_[:,::1], int_[::1])', parallel=True)
def computeTriangle(nodes, neighbourSlices, allNeighbours):
nodeIdToPos = np.empty(np.max(nodes) 1, dtype=np.int_)
for i, node in enumerate(nodes):
nodeIdToPos[node] = i
neighbourSlices, allNeighbours = filterNeighbours(nodes, neighbourSlices, allNeighbours)
threadCount = nb.np.ufunc.parallel.get_num_threads()
nodeSlides = splitSlices(neighbourSlices, threadCount)
s = 0
for threadId in nb.prange(threadCount):
start, stop = nodeSlides[threadId]
for i in range(start, stop):
node1 = nodes[i]
start1, stop2 = neighbourSlices[i]
neighbours1 = allNeighbours[start1:stop2]
for node2 in neighbours1:
start2, stop2 = neighbourSlices[nodeIdToPos[node2]]
neighbours2 = allNeighbours[start2:stop2]
for node3 in neighbours2:
if node3 > node2:
foundId = np.searchsorted(neighbours1, node3)
s = foundId < neighbours1.size and neighbours1[foundId] == node3
return s
# Conversion to Numpy arrays
edgeCount = sum(len(neighbours) for _, neighbours in G.adjacency())
nodes = np.fromiter(G.nodes(), dtype=np.int_)
neighbourSlices = np.empty((len(G), 2), dtype=np.int_)
allNeighbours = np.empty(edgeCount, dtype=np.int_)
curPos = 0
for i, nodeInfos in enumerate(G.adjacency()):
node, edgeInfos = nodeInfos
neighbours = np.fromiter(edgeInfos.keys(), dtype=np.int_)
neighbourSlices[i, 0] = curPos
neighbourSlices[i, 1] = curPos neighbours.size
allNeighbours[curPos:curPos neighbours.size] = neighbours
curPos = neighbours.size
# Actual computation (very fast)
triangleCount = computeTriangle(nodes, neighbourSlices, allNeighbours)
在我的6核機,這種解決方案快45倍比在示例曲線圖的原始迭代解和高達290倍的速度上nx.erdos_renyi_graph(15000, 0.005)。當平均度數(N_edges/N_nodes)很大時,它比其他方法快得多。
不幸的是,大部分時間都花在將圖(順序)轉換為 Numpy 陣列上。因此,如果您想要更快的方法(例如快 2~3 倍),那么您當然需要使用本地語言實作(如 C 或 C )。或者主要使用 Numpy 陣列而不是 NetworkX 圖,但這種方法顯然不太靈活且存檔復雜。
轉載請註明出處,本文鏈接:https://www.uj5u.com/gongcheng/402191.html
上一篇:為什么從生成器運算式生成元組比從串列推導式生成元組慢?
下一篇:理解中的兩個命令
