Finding path with smallest GCD of nodes's weights in directed graph

Question

I wanted to solve this problem: GCD on directed graph

I'm new to SCC, topological ordering, Kosajaru algorithm etc.

Generally, I think that the more nodes we use in path the better result can be, because the cost (GCD) will never grow. This is my approach to this problem:

1. Find all strongly connected components (so I use as many nodes as possible) and the GCD (cost of traversing whole component)
2. We can map found SCCs to nodes and make DAG where single node represents whole SCC

This way we reduced the problem to finding the shortest path (from any to any node) in weighted DGA. But I have no idea how can I approach it with acceptable time complexity. That's how my current code looks:

#include <iostream>
#include <vector>
#include <list>
#include <stack>

using namespace std;

int gcd(int a, int b) {
    if (b == 0)
        return a;

    return gcd(b, a % b);
}

void top_sort(const vector<list<int>>& g, vector<bool>& visited, vector<int>& order, int node) {
    visited[node] = true;

    for (const int ngbr : g[node]) {
        if (visited[ngbr])
            continue;

        top_sort(g, visited, order, ngbr);
    }

    order.push_back(node);
}

void dfs_0(const vector<list<int>>& rg, const vector<int>& c, vector<bool>& visited, vector<int>& components_gcd, vector<int>& components, int node, int component) {
    visited[node] = true;
    components_gcd[component] = gcd(components_gcd[component], c[node]);
    components[node] = component;

    for (const int ngbr : rg[node]) {
        if (visited[ngbr])
            continue;

        dfs_0(rg, c, visited, components_gcd, components, ngbr, component);
    }
}

int solve(int n, vector<int> c, vector<vector<int>> edges) {
    vector<list<int>> g(n + 1, list<int>());
    vector<list<int>> rg(n + 1, list<int>());
    for (int i = 0; i < edges.size(); ++i) {
        g[edges[i][0]].push_back(edges[i][1]);
        rg[edges[i][1]].push_back(edges[i][0]);
    }

    vector<bool> visited(n + 1);
    vector<int> order;
    order.reserve(n);

    for (int i = 1; i <= n; ++i) {
        if (!visited[i])
            top_sort(g, visited, order, i);
    }

    reverse(order.begin(), order.end());
    fill(visited.begin(), visited.end(), false);

    vector<int> components_gcd(n + 1);
    vector<int> components(n + 1);
    int component = 0;
    for (const int node : order) {
        if (!visited[node]) {
            ++component;
            dfs_0(rg, c, visited, components_gcd, components, node, component);
        }
    }

    vector<list<int>> scc_g(component, list<int>());
    for (int i = 0; i < edges.size(); ++i) {
        if (components[edges[i][0]] != components[edges[i][1]]) {
            scc_g[components[edges[i][0]]].push_back(components[edges[i][1]]);
        }
    }

    // I struggle here. scc_g is our DAG with SCCs as single nodes
}

Please comment if something is unclear. Thanks

Answer 1

1. Find all SCCs.
2. Build a DAG.
3. Apply topological sorting.
4. Perform dynamic programming (DP) on nodes in topological order.

Implementation in C++ :

#include <bits/stdc++.h>
using namespace std;

int main(){

    int n,m;
    cin>>n>>m;
    vector<int> cost(n+1,0), adj[n+1], inverted[n+1];
    for(int i=1; i<=n; i++) cin>>cost[i];
    for(int i=0; i<m; i++) {
        int a,b;
        cin>>a>>b;
        adj[a].push_back(b);
        inverted[b].push_back(a);
    }

    vector<int> order, vis(n+1,0);

    function<void(int)> dfs = [&](int node) {
        vis[node] = 1;
        for(auto i: adj[node]) if(!vis[i]) dfs(i);
        order.push_back(node);
    };

    for(int i=1; i<=n; i++) if(!vis[i]) dfs(i);

    reverse(order.begin(), order.end());
    vis.assign(n+1,0);
    
    vector<int> roots(n+1,0);
    vector<vector<int>> components;
    for(int i=0; i<n; i++){
        if(!vis[order[i]]) {
            vector<int> component;
            function<void(int)> dfs2 = [&](int node) {
                vis[node] = 1;
                for(auto i: inverted[node]) if(!vis[i]) dfs2(i);
                component.push_back(node);
            };
            dfs2(order[i]);

            components.push_back(component);
            int root = *min_element(component.begin(), component.end());
            int hcf = cost[order[i]];
            for(auto j: component){
                hcf = __gcd(hcf, cost[j]);
                roots[j] = root;
            }
            cost[root] = hcf;
        }
    }

    vector<int> adjDagInverted[n+1], inDegree(n+1,0), isNode(n+1,0);

    for(int i=1; i<=n; i++) {
        for(auto j: adj[i]) {
            if(roots[i] != roots[j]) {
                isNode[roots[i]] = isNode[roots[j]] = 1;
                adjDagInverted[roots[j]].push_back(roots[i]);
                inDegree[roots[i]]++;  
            }
        }
    }
    
    queue<int> q;
    for(int i=1; i<=n; i++) {
        if(inDegree[i] == 0 && isNode[i]) q.push(i);
    }

    vector<int> topo;
    while(!q.empty()){
        int node = q.front();
        q.pop();
        topo.push_back(node);
        for(auto i: adjDagInverted[node]){
            inDegree[i]--;
            if(inDegree[i] == 0 && isNode[i]) q.push(i);
        }
    }
    
    vector<int> dp(n+1,1e9);
    function<int(int)> dfs3 = [&](int node) -> int{
        if(dp[node] != 1e9) return dp[node];
        int best = cost[node];
        for(auto i: adjDagInverted[node]){
            int nxt = dfs3(i);
            best = min(best, __gcd(cost[node], nxt));
        }
        return dp[node] = best;
    };
    
    int ans = 1e9;
    for(auto it: topo){
        ans = min(ans, dfs3(it)); 
    }
    cout<<ans;
    
    return 0;
}

You can read more about Condensation Graph here : https://cp-algorithms.com/graph/strongly-connected-components.html

Answer 2

Given that:

• GCD(a, b) ≤ a
• GCD(a, b) ≤ b

Then you know that for a graph G(V,E) where the vertex v_n has a cost c_n then for every edge e_(i,j) = (v_i, v_j) the cost c_(i,j) for the edge is GCD(c_i, c_j) which is always less than or equal to both c_i and c_j.

Additionally, GCDs are:

• Reflexive, so GCD(a, b) = GCD(b, a)
• Commutative, so GCD(GCD(a, b), c) = GCD(a, GCD(b, c))

So if you have any graph G(V, E) then:

1. Find the strongly connected components in the graph and:

   • Generate a minor of the graph with each SCC replaced by a single vertex with the cost equal to the GCD of the costs of all the vertices in that SCC.

2. For each vertex with no inbound edges, generate a set of maximal walks starting from that vertex by recursively following all the outbound edges (if they exist) from that vertex to find the maximal walks can calculate the GCD of the cost of all vertices on a path.

   Note: having removed the SCCs there should be no cycles on those paths but you may have to visit vertices multiple times on different paths and branching paths that later converge may have different GCDs. I.e. (12->8->24) and (12->9->24) that share the same start and end have GCDs of 4 and 3 respectively. This means that the algorithm should be O(E) bounded rather than O(V).

3. The GCD for the graph is the minimum of the GCDs in each of those maximal walks.

In python, I think you could implement it using:

from __future__ import annotations

from collections import deque
from io import StringIO
from typing import List, Optional, Tuple, Union


def gcd(a: int, b: int) -> int:
    """Binary GCD algorithm"""
    d: int = 1
    while True:
        a_trailing_zero_bits = (a&-a).bit_length() - 1
        b_trailing_zero_bits = (b&-b).bit_length() - 1
        d <<= min(a_trailing_zero_bits, b_trailing_zero_bits)
        a >>= a_trailing_zero_bits
        b >>= b_trailing_zero_bits
        if a == b:
            return d * a
        elif a < b:
            b = b - a
        else:
            a = a - b


class Graph:
    vertices: List[Vertex]
    edges: List[Edge]
    
    @staticmethod
    def read(stream) -> Graph:
        num_vertices, num_edges = [int(value) for value in stream.readline().split(" ")]
        assert num_vertices > 0
        v = [Vertex(idx, int(cost)) for idx, cost in enumerate(stream.readline().split())]
        assert len(v) == num_vertices
        g = Graph(v)
        for i in range(0, num_edges):
            f, t = [int(index) - 1 for index in stream.readline().split()]
            g.add_edge(Edge(i, v[f], v[t]))
        return g
        
    def __init__(self, vertices: List[Vertex]) -> None:
        self.vertices = vertices
        self.edges = []
       
    def add_edge(self, edge: Edge) -> None:
        self.edges.append(edge)
        
    def __str__(self) -> str:
        vertex_costs = ", ".join(str(v) for v in self.vertices)
        edges = ", ".join(str(e) for e in self.edges)
        return f"G([{vertex_costs}], [{edges}])"
        
    def find_strongly_connected_components(self) -> None:
        """Tarjan's Strongly Connected Component algorithm."""
        index: int = 0
        S: List[Vertex] = []

        def strong_connect(v: Vertex) -> None:
            nonlocal index
            nonlocal S
            v.scc_set_index(index)
            S.append(v)
            index += 1
            
            for e in v.outbound:
                w = e.to_vertex
                if w.scc_index is None:
                    strong_connect(w)
                    assert v.scc_low_index is not None
                    assert w.scc_low_index is not None
                    v.scc_low_index = min(v.scc_low_index, w.scc_low_index)
                elif w.scc_on_stack:
                    assert v.scc_low_index is not None
                    assert w.scc_low_index is not None
                    v.scc_low_index = min(v.scc_low_index, w.scc_low_index)
                
            if v.scc_low_index == v.scc_index:
                w = S.pop()
                w.scc_on_stack = False
                if w is v:
                    return
                scc = StronglyConnectedComponent(w)

                while v is not w:
                    w = S.pop()
                    w.scc_on_stack = False
                    scc.add_vertex(w)

                scc.process_edges()
        
        for v in self.vertices:
            if v.scc_index is None:
                strong_connect(v)

    def find_gcd(self) -> int:
        """Find all Strongly Connected Components and then find all maximal paths
        and return the path calculating the cost as the path is walked and return
        the minimum cost.
        
        The algorithm will terminate early if the minimal cost of 1 is reached.
        """
        assert self.vertices
        min_gcd: Optional[int] = None
        stack: List[Tuple[int, int, Vertex]]
        v: Vertex
        w: Union[Vertex, StronglyConnectedComponent]
        cost: int
        
        self.find_strongly_connected_components()

        for v in self.vertices:
            if v.scc:
                if v.scc.visited or v.scc.inbound:
                    continue
            elif v.inbound:
                continue
            stack = [(v.cost, 0, v)]
            while stack:
                cost, depth, w = stack.pop()
                if w.scc:
                    w.scc.visited = True
                    cost = gcd(cost, w.scc.cost)
                # print(f"{'  '*depth}{w}: {cost}")
                if w.scc:
                    w = w.scc
                if cost == 1:
                    return 1
                elif w.outbound:
                    stack.extend(
                        [
                            (gcd(e.to_vertex.cost, cost), depth + 1, e.to_vertex)
                            for e in w.outbound
                        ],
                    )
                elif min_gcd is None:
                    min_gcd = cost
                else:
                    min_gcd = min(min_gcd, cost)
        assert min_gcd is not None
        return min_gcd
                   
            
   
class Vertex:
    scc_index: Optional[int]
    scc_low_index: Optional[int]
    scc_on_stack: bool
    scc: Optional[StronglyConnectedComponent]
    index: int
    _cost: int
    inbound: List[Edge]
    outbound: List[Edge]
    
    def __init__(self, index: int, cost: int) -> None:
        self.index = index
        self._cost = cost
        self.scc_index = None
        self.scc_low_index = None
        self.scc_on_stack = False
        self.scc = None
        self.inbound = []
        self.outbound = []

    def scc_set_index(self, index: int) -> None:
        self.scc_index = index
        self.scc_low_index = index
        self.scc_on_stack = True

    @property    
    def cost(self) -> int:
        return self.scc.cost if self.scc else self._cost

    def __str__(self) -> str:
        if self.scc:
            return f"{self._cost}({self.scc.cost})"
        else:
            return f"{self._cost}"


class Edge:
    index: int
    from_vertex: Vertex
    to_vertex: Vertex
    scc: Optional[StronglyConnectedComponent]

    def __init__(self, index: int, f: Vertex, t: Vertex) -> None:
        self.index = index
        self.from_vertex = f
        self.to_vertex = t
        self.from_vertex.outbound.append(self)
        self.to_vertex.inbound.append(self)
        self.scc = None
        
    def __str__(self) -> str:
        return f"({self.from_vertex.index} -> {self.to_vertex.index})"
        
class StronglyConnectedComponent:
    vertices: List[Vertex]
    inbound: List[Edge]
    outbound: List[Edge]
    internal: List[Edge]
    visited: bool
    cost: int
    
    def __init__(self, vertex: Vertex) -> None:
        self.vertices = []
        self.inbound = []
        self.outbound = []
        self.internal = []
        self.visited = False
        self.add_vertex(vertex)
        
        
    def add_vertex(self, vertex: Vertex) -> None:
        if self.vertices:
            self.cost = gcd(self.cost, vertex.cost)
        else:
            self.cost = vertex.cost
        self.vertices.append(vertex)
        vertex.scc = self

    def process_edges(self) -> None:
        for v in self.vertices:
            for e in v.inbound:
                 if e.from_vertex.scc is self:
                     if e.scc is None:
                         self.internal.append(e)
                         e.scc = self
                 else:
                     self.inbound.append(e)

            for e in v.outbound:
                 if e.to_vertex.scc is self:
                     if e.scc is None:
                         self.internal.append(e)
                         e.scc = self
                 else:
                     self.outbound.append(e)

for input, expected_gcd in (
    (
        "7 6\n4 6 8 20 12 6 5\n1 2\n2 3\n3 1\n2 4\n3 6\n6 5",
        2,
    ),
    (
        "3 2\n12 16 5\n1 2\n2 1",
        4,
    ),
    (
        "3 2\n8 16 5\n1 2\n2 1",
        5,
    ),
    (
        "3 3\n8 16 5\n1 2\n2 1\n2 3",
        1,
    ),
    (
        "7 8\n24 16 12 40 10 15 100\n1 2\n1 3\n2 4\n3 4\n4 5\n4 6\n5 7\n6 7",
        1,
    ),
    (
        "7 8\n24 16 12 40 10 30 100\n1 2\n1 3\n2 4\n3 4\n4 5\n4 6\n5 7\n6 7",
        2,
    ),
    (
        "7 8\n24 16 12 40 20 60 100\n1 2\n1 3\n2 4\n3 4\n4 5\n4 6\n5 7\n6 7",
        4,
    ),
    (
        # All vertices have inbound edges.
        "6 7\n18 36 24 36 54 27\n1 2\n2 3\n3 1\n4 5\n5 6\n6 4\n2 5",
        3,
    ),
    (
        # No edges.
        "4 0\n5 4 2 3",
        2,
    ),
):
    graph = Graph.read(StringIO(input))
    min_gcd = graph.find_gcd()
    assert min_gcd == expected_gcd
    print(graph, min_gcd)