Multi-Source BFS, in depth: seeding the queue, nearest-source labels, closest pairs of sources and reversed graphs

One wave from many starting points

Add a virtual vertex S to the graph and connect it to every source with an edge of cost zero. The distance from S to any vertex v is then the minimum over sources s of the distance from s to v, which is exactly the quantity we want. A BFS from S would pop S, discover every source at distance zero, and continue. Seeding the queue with all sources is that same search with the first step already taken, so S never needs to exist in memory.

The correctness of plain BFS rests on one invariant: the queue holds vertices in non-decreasing order of distance, and at any moment it contains at most two distinct distance values, d and d + 1. A seeded queue starts with every entry at distance 0, so the invariant holds from the first pop. Every vertex is therefore assigned its distance the first time it is discovered, and that value is final. The total cost is O(V + E), independent of how many sources there are.

Compare that with the alternatives. k separate searches cost O(k(V + E)), and taking a minimum over k distance arrays also needs O(kV) memory unless you merge as you go. Dijkstra from a super-source works but pays a log factor for a heap you do not need when every edge has the same weight. Seeding the queue is the cheapest correct option whenever edges are unweighted.

The template: distance, label and parent

The template carries three arrays: the distance, the source that reached each vertex first, and the parent on that shortest path. The label is what turns a distance map into an assignment, and the parent lets you print the route.

from collections import deque

def multi_source_bfs(adj, sources):
    """adj: list of neighbour lists, vertices 0..n-1. sources: vertex ids.
    Returns dist (hops to nearest source, -1 if unreachable),
    label (the source that reached each vertex) and parent."""
    n = len(adj)
    dist = [-1] * n
    label = [-1] * n
    parent = [-1] * n
    q = deque()
    for s in sources:
        if dist[s] == -1:            # duplicate sources are harmless
            dist[s] = 0
            label[s] = s
            q.append(s)
    while q:
        u = q.popleft()
        for v in adj[u]:
            if dist[v] == -1:        # mark when enqueued, never when popped
                dist[v] = dist[u] + 1
                label[v] = label[u]
                parent[v] = u
                q.append(v)
    return dist, label, parent

def route_to_source(v, parent):
    path = [v]
    while parent[path[-1]] != -1:
        path.append(parent[path[-1]])
    return path                      # v ... source

Two details matter. Marking a vertex when it is enqueued, rather than when it is popped, keeps each vertex in the queue at most once; marking on pop lets the same vertex enter many times and, with many sources, the queue can grow far beyond V. And the duplicate check on the seeds matters when the source list comes from data, where the same depot can appear twice.

Worked example: nearest depot on a road graph

Take nine junctions laid out as a three by three block, A to I, with the eleven roads shown in the diagram (the road between B and E is closed). Depots sit at A and I. The seeds go in as [A, I]. The trace, one pop per row:

Every vertex was visited once and nine pops answered the question for both depots. Notice C, E and G. Each is exactly two hops from both depots, and all three were labelled A. Nothing about the graph made A nearer; A simply sat ahead of I in the queue, so its wave reached the tied vertices first. Reverse the seed order and all three flip to I. If a downstream system depends on those labels, for example to split delivery work, that order dependence is a bug waiting to be reported.

Nearest-source regions and ties

The labels partition the reachable graph into regions, the graph analogue of a Voronoi diagram. Distances are unique; labels at ties are not. You have three honest options.

• Accept arbitrary ties and document that they follow seed order. Fine for a game AI choosing which enemy to chase.
• Break ties by a rule, such as the smallest source id. Process the search level by level: build the whole next frontier, and when a vertex at distance d + 1 is discovered again from another vertex at distance d, keep the smaller label. Because all level-d labels are final before level d + 1 is built, the rule propagates correctly.
• Keep every tied source by storing a small set per vertex. This costs memory, but it is the right answer when a tie means shared responsibility, such as two clinics equally close to a village.

def multi_source_bfs_min_label(adj, sources):
    n = len(adj)
    dist, label = [-1] * n, [-1] * n
    frontier = sorted(set(sources))
    for s in frontier:
        dist[s], label[s] = 0, s
    d = 0
    while frontier:
        nxt = []
        for u in frontier:
            for v in adj[u]:
                if dist[v] == -1:
                    dist[v], label[v] = d + 1, label[u]
                    nxt.append(v)
                elif dist[v] == d + 1 and label[u] < label[v]:
                    label[v] = label[u]      # same level, smaller source wins
        frontier, d = nxt, d + 1
    return dist, label

The closest pair of sources from one pass

A question that looks harder is also answered by the labels: which two sources are closest to each other, and how far apart are they? Running BFS from every source costs O(k(V + E)). After one multi-source pass, scan every edge (u, v) whose endpoints carry different labels and take the minimum of dist[u] + 1 + dist[v].

def closest_source_pair(adj, dist, label):
    best, pair = float("inf"), None
    for u in range(len(adj)):
        if dist[u] == -1:
            continue
        for v in adj[u]:
            if label[v] != -1 and label[u] != label[v]:
                d = dist[u] + 1 + dist[v]
                if d < best:
                    best, pair = d, (label[u], label[v])
    return best, pair

Why it is exact: every candidate is the length of a real walk from label[u] to u, across the edge, and on to label[v], which joins two different sources, so no candidate is shorter than the true closest distance L. In the other direction, take a shortest path of length L between the closest pair a and b. Its first vertex is labelled a and its last b, so some edge (u, v) along it changes label. dist[u] is at most the distance from a to u along the path, and dist[v] is at most the distance from v to b, so that edge's candidate is at most L. The minimum is therefore exactly L. In the example the four red edges all give 2 + 1 + 1 = 4, matching the route A, B, C, F, I.

Two narrower questions fall out of the same pass. For the distance from set X to set Y, seed X and read the smallest dist over Y; you can stop the search the first time a Y vertex is discovered. For weighted graphs, replace the + 1 with the edge weight and the BFS with a seeded Dijkstra; the proof is unchanged.

Directed graphs: search from the targets

On a directed graph, BFS from the sources computes the distance from the nearest source to each vertex. Many real questions run the other way: from every service, how many calls until it reaches any database? From every web page, how many clicks to any checkout page? Those are distances to the nearest target, and you get them by seeding the targets and searching the reversed graph, following each edge backwards. Building the reverse adjacency list is one O(V + E) pass. Forgetting it gives answers that look plausible and are wrong for every vertex whose in-edges and out-edges differ.

The same reversal underlies the boundary-flood family: to find which cells can drain to an edge, seed the edge cells and expand to neighbours that could flow into them. One search from the destination replaces one search per start point.

Staggered start times and changing sources

Sometimes sources do not all start at time zero: a fire breaks out at A at minute 0 and at I at minute 3, or warehouses open at different hours. Seeding everything at once is now wrong, because I would spread three minutes too early. The fix is to keep the queue ordered. With small integer offsets, sort the seeds by start time and inject each one when the frontier reaches its time, before expanding that level; a late source that an earlier wave already reached is simply skipped. With large or real-valued offsets, use Dijkstra from the super-source with edge weights equal to the offsets. If edges also have weights 0 and 1, the same idea works with a deque, as in 0-1 BFS.

Adding a source later is cheap. Run a BFS from the new source alone, but only push a vertex when the new distance is strictly smaller than the stored one; the search stops where the old regions are already closer. Removing a source is not cheap: vertices in its region may now belong to several neighbours, and the simplest correct approach is to recompute, or at least to clear its region and re-seed from that region's boundary with the distances of the surviving neighbours, injected in distance order exactly as for staggered starts.

Failure modes

The bugs that recur in code review:

• Marking on pop. Vertices enter the queue many times; correctness survives but memory and time can explode on dense graphs.
• Seeding one source at a time in a loop with a fresh visited array. This is k searches, not one, and it silently turns O(V + E) into O(k(V + E)).
• Forgetting unreachable vertices. A -1 distance must be checked before it is used in arithmetic; the closest-pair scan above skips them for that reason.
• Searching the wrong direction on a directed graph, covered above.
• Relying on tie labels that depend on seed order, including the order a database returned the sources in.
• Staggered starts seeded at once, which makes late sources too strong.

Running it on large graphs

At millions of vertices the algorithm is fine and the data structures are the problem. Store the graph in compressed sparse row form: one offsets array of length V + 1 and one targets array of length E. Use int32 arrays for dist, label and parent; at 100 million vertices each costs 400 MB, so drop parent unless you need routes, and drop label unless you need the assignment. Replace the deque with two flat frontier arrays, which is also the shape parallel BFS needs. On low-diameter graphs, direction-optimising BFS, where late levels scan unvisited vertices for a visited neighbour instead of expanding the frontier, cuts edge inspections sharply, and multiple sources make the frontier large early, which is when that switch pays.

For the general BFS foundations, including edge classification and the visited-set argument, see BFS and DFS in depth. For searches whose states are not plain vertices, such as a position plus remaining fuel, see BFS on implicit graphs; multi-source seeding works the same way there, by seeding every start state.

What to do next

1. Write the multi_source_bfs function above and test it against k single-source searches plus a minimum on random graphs; the two must agree on every distance.
2. Decide what a tie means in your problem and pick one of the three tie policies explicitly.
3. On a directed graph, write down whether you need distance from or distance to the sources, and build the reverse adjacency list if it is the latter.
4. If you need the two closest sources, add the edge scan and check it on a small graph by brute force.
5. If sources start at different times, sort and inject them by level, or switch to Dijkstra from a virtual source.
6. For graphs over ten million vertices, move to CSR arrays and int32 state, and measure memory before speed.