CommonFlowDecomp.py

import networkx as nx
import kCommonFlowDecomp as kCFD
from typing import Optional
from itertools import count
import time
import copy
import utils
from collections import deque
from queue import Queue
from typing import Optional


class AbstractSourceSinkGraph(nx.DiGraph):
    """Base class for s-t augmented graphs (internal).

    This class (introduced when unifying :class:`stDAG` and :class:`stDiGraph`) factors out
    logic previously duplicated in both classes.

    Core responsibilities
    ---------------------
    * Store and expose the original ``base_graph`` plus user supplied ``additional_starts`` / ``additional_ends``.
    * Validate that all nodes are strings and that any additional start/end nodes belong to ``base_graph``.
    * Create unique global source / sink node identifiers (``self.source`` / ``self.sink``).
    * Attach the global source to every (in-degree 0) source or additional start; attach every (out-degree 0) sink
      or additional end to the global sink.
    * Expose convenience collections: ``source_edges``, ``sink_edges``, ``source_sink_edges``.
    * Provide shared flow helper utilities: :meth:`get_non_zero_flow_edges` and
      :meth:`get_max_flow_value_and_check_non_negative_flow`.

    Extension hooks
    ---------------
    Subclasses customise behaviour via two lightweight hooks:

    * ``_pre_build_validate`` - extra validation before augmentation (e.g. acyclicity for ``stDAG``).
    * ``_post_build`` - populate subclass specific derived structures (e.g. condensation for ``stDiGraph``).

    Backwards compatibility
    -----------------------
    External code should keep instantiating [`stDAG`](stdag.md) or [`stDiGraph`](stdigraph.md); their public APIs
    are unchanged. ``AbstractSourceSinkGraph`` is an internal implementation detail and may change without notice.
    """

    def __init__(
        self,
        base_graph: nx.DiGraph,
        additional_starts: Optional[list] = None,
        additional_ends: Optional[list] = None,
    ):
        if not all(isinstance(node, str) for node in base_graph.nodes()):
            # utils.logger.error(f"{__name__}: Every node of the graph must be a string.")
            raise ValueError("Every node of the graph must be a string.")

        super().__init__()
        self.base_graph = base_graph
        if "id" in base_graph.graph:
            self.id = str(base_graph.graph["id"])
        else:
            self.id = str(id(self))

        self.additional_starts = set(additional_starts or [])
        self.additional_ends = set(additional_ends or [])

        # Ensure any declared additional start/end nodes are in the base graph
        if not self.additional_starts.issubset(base_graph.nodes()):
            # utils.logger.error(f"{__name__}: Some nodes in additional_starts are not in the base graph.")
            raise ValueError("Some nodes in additional_starts are not in the base graph.")
        if not self.additional_ends.issubset(base_graph.nodes()):
            # utils.logger.error(f"{__name__}: Some nodes in additional_ends are not in the base graph.")
            raise ValueError("Some nodes in additional_ends are not in the base graph.")

        self.source = f"source_{id(self)}"
        self.sink = f"sink_{id(self)}"

        # Hooks
        self._pre_build_validate()
        self._augment_with_source_sink()
        self._post_build()

        nx.freeze(self)

    # ----------------------------- Hooks ---------------------------------
    def _pre_build_validate(self):  # pragma: no cover - default is no-op
        pass

    def _post_build(self):  # pragma: no cover - default is no-op
        pass

    # --------------------------- Build logic -----------------------------
    def _augment_with_source_sink(self):
        # Add base nodes/edges
        self.add_nodes_from(self.base_graph.nodes(data=True))
        self.add_edges_from(self.base_graph.edges(data=True))

        # Connect global source & sink
        for u in self.base_graph.nodes:
            if self.base_graph.in_degree(u) == 0 or u in self.additional_starts:
                self.add_edge(self.source, u)
            if self.base_graph.out_degree(u) == 0 or u in self.additional_ends:
                self.add_edge(u, self.sink)

        self.source_edges = list(self.out_edges(self.source))
        self.sink_edges = list(self.in_edges(self.sink))
        self.source_sink_edges = set(self.source_edges + self.sink_edges)

    # ----------------------- Shared helper methods -----------------------
    def get_non_zero_flow_edges(
        self, flow_attr: str, edges_to_ignore: set = set()
    ) -> set:
        """Return set of edges whose attribute `flow_attr` is non-zero and not ignored."""
        non_zero_flow_edges = set()
        for u, v, data in self.edges(data=True):
            if (u, v) not in edges_to_ignore and data.get(flow_attr, 0) != 0:
                non_zero_flow_edges.add((u, v))
        return non_zero_flow_edges

    def get_max_flow_value_and_check_non_negative_flow(
        self, flow_attr: str, edges_to_ignore: set
    ) -> float:
        """Return maximum value of `flow_attr` over edges (ignoring some) verifying non-negativity.

        Raises ValueError if any required attribute missing or negative.
        """
        w_max = float("-inf")
        if edges_to_ignore is None:
            edges_to_ignore = set()
        for u, v, data in self.edges(data=True):
            if (u, v) in edges_to_ignore:
                continue
            if flow_attr not in data:
                # utils.logger.error(f"Edge ({u},{v}) does not have the required flow attribute '{flow_attr}'.")
                raise ValueError(
                    f"Edge ({u},{v}) does not have the required flow attribute '{flow_attr}'."
                )
            if data[flow_attr] < 0:
                # utils.logger.error(f"Edge ({u},{v}) has negative flow value {data[flow_attr]}. All flow values must be >=0.")
                raise ValueError(
                    f"Edge ({u},{v}) has negative flow value {data[flow_attr]}. All flow values must be >=0."
                )
            w_max = max(w_max, data[flow_attr])
        return w_max

class Arc_Dominator_Tree:

    def __init__(self, n:int, start:str, idoms:dict, edgelist : list, X : set, id:str):
        self.id              = id
        self.n               = n
        self.start           = start
        self.X               = X
        self.idom            = idoms
        self.children        = {e: [] for e in edgelist}
        self.children[start] = []

        for node,idom in self.idom.items():
            self.children[idom].append(node)

        self.idom_X            = dict()
        # self.idom_X[start]     = start
        self.children_X        = {e: [] for e in X}
        self.children_X[start] = []
        self.build_children_relation_X()

    def is_leaf_X(self, arc : tuple):
        # if arc not in self.children_X:
        #     return False
        return len(self.children_X[arc])==0

    def has_unique_child_X(self, arc : tuple):
        # if arc not in self.children_X:
        #     return False
        return len(self.children_X[arc])==1
    
    def get_dominators(self, arc : tuple):
        dominators = []
        while arc != self.start:
            dominators.append(arc)
            arc = self.idom[arc]
        return dominators
    
    def build_children_relation_X(self):

        def dfs(node, last_in_X): # recall that X is a set of arcs. the term "node" is to allude to nodes of the dominator tree
            if node != last_in_X and node in self.X: # note that sink and source are never in X
                self.children_X[last_in_X].append(node)
                self.idom_X[node] = last_in_X
                last_in_X = node
            for child in self.children[node]:
                dfs(child, last_in_X)

        dfs(self.start, self.start)

    #a unitary path in a dominator tree is a path towards the root such that every node has exactly one children except the deepest node
    def find_unitary_path_X(self, arc : tuple, mode : str):
        if mode == "up":
            fn = ( lambda node : self.idom_X[node]        if self.has_unique_child_X(self.idom_X[node]) and self.idom_X[node] != self.start else node )
        if mode == "down":
            fn = ( lambda node : self.children_X[node][0] if self.has_unique_child_X(node)                                                  else node )

        path = [arc]
        while arc != fn(arc):
            arc = fn(arc)
            path.append(arc)
        return path

class CommonFlowDecomp:
    def __init__(self, G: nx.DiGraph, num_flows: int, flow_attr: str = "flow", subpath_constr=None, console=False, greedy_width_solve=False, use_safe_sequences=False, use_safe_paths=False):        
        self.G = G
        
        self.num_flows = num_flows
        # Can we use a better lower bound?
        self.maximum_k = self.G.number_of_edges() * num_flows
        self.flow_attr = flow_attr
        self.console = console
        self.source = next(v for v in self.G.nodes() if self.G.in_degree(v) == 0)
        self.sink = next(v for v in self.G.nodes() if self.G.out_degree(v) == 0)

        self.subpath_constr = subpath_constr or []

        self.greedy_width_solve = greedy_width_solve
        
        # Used for optimizations
        self.safe_sequences = []
        
        # Logging Data
        self.safe_sequence_time = None
        self.safe_paths_time = None
        self.safe_paths_store = []
        self.safe_sequences_store = []
        
        if use_safe_sequences:
            t0 = time.perf_counter()
            aug_G = AbstractSourceSinkGraph(self.G)
            safe_seqs = self.maximal_safe_sequences_via_dominators(aug_G, set(self.G.edges()))
            self.safe_sequences_store = self.clean_safe_sequences(safe_seqs)
            self.safe_sequences += self.safe_sequences_store
            self.safe_sequence_time = time.perf_counter() - t0
        
        if use_safe_paths:
            t0 = time.perf_counter()
            def take_path_single(flow_network, path, weight, i):
                for j in range(len(path) - 1):
                    flow_network[path[j]][path[j+1]]['flow'][i] -= weight

            def is_decomposed_single(flow_network, i, source):
                for u, v in flow_network.edges(source):
                    if flow_network[u][v]['flow'][i] != 0:
                        return False
                return True

            def find_path_with_highest_flow_single(graph, start, end, i):
                vertices = list(nx.topological_sort(graph))
                best_flow = {v: -1 for v in vertices}
                best_path = {v: None for v in vertices}
                
                best_flow[start] = float("inf")
                best_path[start] = (start,)
                
                for v in vertices:
                    if v == start:
                        continue
                    for u in graph.predecessors(v):
                        flow = graph[u][v]["flow"][i]
                        incoming = min(best_flow[u], flow)
                        if incoming > best_flow[v]:
                            best_flow[v] = incoming
                            best_path[v] = best_path[u] + (v,)
                
                return list(best_path[end]), best_flow[end]

            def greedily_decompose_flow_single(flow_network, i, source):
                paths = []
                weights = []
                sink = next(v for v in flow_network.nodes() if flow_network.out_degree(v) == 0)
                while not is_decomposed_single(flow_network, i, source):
                    path, weight = find_path_with_highest_flow_single(flow_network, source, sink, i)
                    take_path_single(flow_network, path, weight, i)
                    paths.append(path)
                    weights.append(weight)
                return paths, weights
        
            decomp_paths = []
            
            g_copy = copy.deepcopy(self.G)
            for i in range(num_flows):
                paths, _ = greedily_decompose_flow_single(g_copy, i, '0')
                decomp_paths += paths

            self.safe_paths_store = self.compute_union_flow_safe_paths(G=self.G, flow_attr="flow", decomp_paths=decomp_paths, num_flows=num_flows, no_duplicates=True)
            
            self.safe_sequences += self.safe_paths_store
            self.safe_paths_time = time.perf_counter() - t0
        
        
    def compute_union_flow_safe_paths(self, G: nx.DiGraph, flow_attr: str, decomp_paths, num_flows, no_duplicates: bool = True):
        unique_safe_paths = []
        unique_paths = set()
            
        for i in range(num_flows):
            safe_paths = self.compute_flow_decomp_safe_paths(G=G, flow_attr=flow_attr, decomp_paths=decomp_paths, flow_index=i, no_duplicates=no_duplicates)

            for p in safe_paths:
                if not tuple(p) in unique_paths:
                    unique_paths.add(tuple(p))
                    unique_safe_paths.append(p)
                    
        return unique_safe_paths

    def compute_inexact_flow_decomp_safe_paths(self, 
        G: nx.DiGraph, 
        lowerbound_attr: str, 
        upperbound_attr: str, 
        decomp_paths: list, 
        flow_index,
        no_duplicates: bool = True
    ) -> list:
        """
        Computes all flow-decomposition safe paths for a given non-negative inexact flow,
        given as intervals [lowerbound_attr, upperbound_attr] for every edge.
        A path is called *flow-decomposition safe* if for all flow decompositions of an inexact flow,
        it appears as a subpath of some path of the decomposition.
        See https://doi.org/10.1109/BIBM47256.2019.8983180, https://doi.org/10.1007/978-3-031-04749-7_14, 
        https://doi.org/10.4230/LIPIcs.SEA.2024.14

        Parameters
        ----------
        - `G`: nx.DiGraph: 

            A directed graph as [networkx DiGraph](https://networkx.org/documentation/stable/reference/classes/digraph.html).

        - `lowerbound_attr`: str

            The name of the edge attribute where to get the lower bound flow values from.

        - `upperbound_attr`: str

            The name of the edge attribute where to get the upper bound flow values from.

        - `decomp_paths`: list

            The list of paths from which flow decomposition safe paths will be extracted.
            It is recommended that the paths correspond to a flow decomposition, otherwise the output
            does not necessarily correspond to flow decomposition safe paths.

        - `no_duplicates`: bool

            If `True`, the function returns a set of paths without duplicates.
            
        Returns
        -------
        
        - `paths` (list of lists):
        
            A list of maximal-length flow decomposition safe paths, as lists of edges.

        Raises
        ------
        - ValueError

            - If an edge in a path from decomp_paths does not have the required flow attributes `lowerbound_attr` and `upperbound_attr`.
            - If an edge in a path from decomp_paths has a negative flow lower or upper bound.
            - If an edge in a path from decomp_paths has a larger lower bound than upper bound flow attribute.
            - If an edge in a path from decomp_paths has a zero flow upper bound.
        """
        # Check the necessary constraints
        for path in decomp_paths:
            for u, v in zip(path, path[1:]):
                for flow_attr in [lowerbound_attr, upperbound_attr]:
                    if flow_attr not in G.edges[u, v]:
                        utils.logger.error(
                            f"{__name__}: Edge ({u},{v}) does not have the required flow attribute '{flow_attr}'. Check that the attribute passed under 'flow_attr' is present in the edge data."
                        )
                        raise ValueError(
                            f"Edge ({u},{v}) does not have the required flow attribute '{flow_attr}'. Check that the attribute passed under 'flow_attr' is present in the edge data."
                        )
                if G.edges[u, v][lowerbound_attr][flow_index] < 0:
                    raise ValueError(
                        f"Edge ({u},{v}) has negative lower bound flow value {G.edges[u, v][lowerbound_attr][flow_index]}. All lower bound flow values must be >=0."
                    )
                if G.edges[u, v][lowerbound_attr][flow_index] > G.edges[u, v][upperbound_attr][flow_index]:
                    raise ValueError(
                        f"Edge ({u},{v}) has a larger lower bound flow value {G.edges[u, v][lowerbound_attr][flow_index]} than upper bound flow value {G.edges[u, v][upperbound_attr][flow_index]}."
                    )


        safe_paths_set = set()
        safe_paths_list = []

        # The algorithm follows a two pointer approach computing inexact excess flow
        # See https://doi.org/10.1007/978-3-031-04749-7_11 and https://doi.org/10.4230/LIPIcs.SEA.2024.14

        for path in decomp_paths:
            if len(path) <= 1:
                continue

            safe_path = deque()
            L, R = 0, 0
            inexact_excess = 0
            safe_path.append(path[L])
            path_not_suffix_of_previous = True

            while R+1 < len(path):
                # Initialize new safe path
                if L == R:
                    assert len(safe_path) == 1
                    assert inexact_excess == 0

                    R += 1
                    inexact_excess = G.edges[path[L], path[R]][lowerbound_attr][flow_index]
                    safe_path.append(path[R])
                    path_not_suffix_of_previous = True

                # Maximally extend the safe path to the right
                while R+1 < len(path):
                    rightdiff = G.edges[path[R], path[R+1]][upperbound_attr][flow_index] - sum(G.edges[u, v][upperbound_attr][flow_index] for u, v in G.out_edges(path[R]))

                    if inexact_excess + rightdiff <= 0:
                        break

                    inexact_excess += rightdiff
                    safe_path.append(path[R+1])
                    R += 1
                    path_not_suffix_of_previous = True

                if path_not_suffix_of_previous:
                    safe_paths_set.add(tuple(safe_path.copy())) if no_duplicates else safe_paths_list.append(safe_path.copy())

                # Remove the left most edge of the safe path
                inexact_excess -= G.edges[path[L], path[L+1]][lowerbound_attr][flow_index]
                if L+1 < R:
                    inexact_excess += sum(G.edges[u, v][upperbound_attr][flow_index] for u, v in G.out_edges(path[L+1])) - G.edges[path[L+1], path[L+2]][upperbound_attr][flow_index]
                    inexact_excess += G.edges[path[L+1], path[L+2]][lowerbound_attr][flow_index]
                safe_path.popleft()
                L += 1
                path_not_suffix_of_previous = False

        if no_duplicates:
            safe_paths_list = [list(sp) for sp in safe_paths_set]
        
        # converting each path from a list of nodes to a list of edges
        safe_paths_list_edges = []
        for safe_path in safe_paths_list:
            assert len(safe_path) > 1
            safe_path_edges = []
            for i in range(len(safe_path)-1):
                safe_path_edges.append((safe_path[i], safe_path[i+1]))
            safe_paths_list_edges.append(safe_path_edges)

        return safe_paths_list_edges

    def compute_flow_decomp_safe_paths(self, 
        G: nx.DiGraph, 
        flow_attr: str, 
        decomp_paths,
        flow_index,
        no_duplicates: bool = True
    ) -> list:
        """
        Computes all flow-decomposition safe paths for a given non-negative flow.
        A path is called *flow-decomposition safe* if for all flow decompositions,
        it appears as a subpath of some path of the decomposition.

        Parameters
        ----------
        - `G`: nx.DiGraph:

            A directed graph as [networkx DiGraph](https://networkx.org/documentation/stable/reference/classes/digraph.html).

        - `flow_attr`: str

            The name of the edge attribute where to get the flow values from.

        - `no_duplicates`: bool

            If `True`, the function returns a set of paths without duplicates.

        Returns
        -------

        - `paths` (list of lists):

            A list of flow safe paths, as lists of edges.

        Raises
        ------

        - ValueError: If an edge does not have the required flow attribute.
        - ValueError: If an edge has a negative flow value.
        """

        # stG = stdag.stDAG(G)
        # decompose it here
        # decompose each of the flows
        # decomp_paths = stG.decompose_using_max_bottleneck(flow_attr)[0]
        return self.compute_inexact_flow_decomp_safe_paths(
            G = G, 
            lowerbound_attr = flow_attr, 
            upperbound_attr = flow_attr, 
            decomp_paths = decomp_paths, 
            flow_index=flow_index,
            no_duplicates = no_duplicates,
            )

    def find_idom(self, adj_dict, s, t) -> list:
        
        p = self.find_path(adj_dict, s, t)


        for i in range(len(p)-1):
            adj_dict[p[i]].remove(p[i+1])
            adj_dict[p[i+1]].append(p[i])

        n            = len(adj_dict)
        i            = 1
        component = dict()  # [0] * n
        for v in adj_dict.keys():
            component[v] = 0
        q            = Queue(maxsize = n+1)
        component[s] = 1
        first_node   = 0
        first_bridge = None
        q.put(s)

        while component[t]==0:

            if i!=1:
                while component[p[first_node]] != 0:
                    first_node += 1

                first_bridge = ( p[first_node-1] ,p[first_node] )
                break

            while not q.empty():
                u = q.get()
                for v in adj_dict[u]:
                    if component[v]==0:
                        q.put(v)
                        component[v]=i
            i = i+1

        for i in range(len(p)-1):
            u,v = p[i],p[i+1]
            adj_dict[v].pop()
            adj_dict[u].append(v)

        return first_bridge

    def find_path(self, adj_dict, s, t):
        """Find a path from s to t using DFS."""
        def dfs_path(node, path: list, visited: set):
            if node == t:
                return True
            visited.add(node)
            for neighbor in adj_dict[node]:
                if neighbor not in visited:
                    path.append(neighbor)
                    if dfs_path(neighbor, path, visited):
                        return True
                    path.pop()  # Backtrack if this path doesn't lead to t
            return False
        
        path = [s]
        visited = set()
        dfs_path(s, path, visited)
        return path

    def maximal_safe_sequences_via_dominators(self, G : AbstractSourceSinkGraph, X = set()) -> list :
        if X == None or len(X) == 0:
            return []

        s_idoms = dict()
        t_idoms = dict()

        adj_dict = {u: list(G.successors(u)) for u in G.nodes()}
        adj_dict_rev = {u: list(G.predecessors(u)) for u in G.nodes()}

        for (u,v) in G.edges:

            s_idom = self.find_idom(adj_dict_rev, u, G.source)
            t_idom = self.find_idom(adj_dict    , v,   G.sink)
            s_idoms[(u,v)] = tuple(reversed(s_idom)) if s_idom != None else G.source
            t_idoms[(u,v)] = t_idom                  if t_idom != None else G.sink

        T_s = Arc_Dominator_Tree(G.number_of_nodes(), G.source, s_idoms, G.edges, X, G.id+str("_s-domtree"))
        T_t = Arc_Dominator_Tree(G.number_of_nodes(), G.sink  , t_idoms, G.edges, X, G.id+str("_t-domtree"))

        leaves_s_X = [ node for node,children in T_s.children_X.items() if len(children)==0 ]
        cores = []
        for leaf in leaves_s_X:
            s_unitary_path = T_s.find_unitary_path_X(leaf,   "up")
            t_unitary_path = T_t.find_unitary_path_X(leaf, "down")
            
            if len(t_unitary_path) > len(s_unitary_path):
                continue

            i=0
            good_sequence = True
            while good_sequence and i<len(t_unitary_path):
                if s_unitary_path[i] != t_unitary_path[i]:
                    good_sequence = False
                else:
                    i=i+1
                    
            if i == len(t_unitary_path) and T_t.is_leaf_X(t_unitary_path[len(t_unitary_path)-1]):
                assert(good_sequence==True)
                cores.append(leaf)

        maximal_safe_sequences = []
        for core in cores: # O(length of all maximal safe sequences), with no duplicates
            s_doms = T_s.get_dominators(core)
            t_doms = T_t.get_dominators(core)
            maximal_safe_sequences.append(s_doms[::-1] + t_doms[1:] )
        
        return maximal_safe_sequences

    def clean_safe_sequences(self, safe_sequences):
        new_sequences = []
        for s in safe_sequences:
            new_sequences.append(s[1:-1])
        return new_sequences
     
    def maximal(self, profiles):
        dominant = set()
        for p, path_p in profiles:
            dominated = False
            for q, _ in profiles:
                if all(q[i] >= p[i] for i in range(len(p))) and any(q[i] > p[i] for i in range(len(p))):
                    dominated = True
                    break
            if not dominated:
                dominant.add((p, path_p))
        return dominant

    def find_path_with_highest_flow(self):
        vertices = list(nx.topological_sort(self.G))
        bottle_necks = {v: set() for v in vertices}
        bottle_necks[self.source] = {(tuple([float("inf")] * self.num_flows), (self.source,))}
        for v in vertices:
            if v == self.source:
                continue
            s = set()
            for u in self.G.predecessors(v):
                flows = self.G[u][v]["flow"]
                for profile, path in bottle_necks[u]:
                    new_profile = tuple(min(profile[i], flows[i]) for i in range(self.num_flows))
                    s.add((new_profile, path + (v,)))
            bottle_necks[v] = self.maximal(s)
        # end = str(max(bottle_necks.keys(), key=int))
        best_profile, best_path = max(bottle_necks[self.sink], key=lambda x: sum(x[0]))
        return list(best_path), list(best_profile)

    def take_path(self, path, weights):
        for i in range(0, len(path)-1):
            for j in range(0, len(weights)):
                self.G[path[i]][path[i+1]]['flow'][j] -= weights[j]

    def greedily_decompose_flow(self):
        paths = []
        weights = []
        sink = next(v for v in self.G.nodes() if self.G.out_degree(v) == 0)
        i = 0
        while not utils.is_decomposed(self.G, self.num_flows, self.source):
            path, weight = self.find_path_with_highest_flow()
            self.take_path(path, weight)
            paths.append(path)
            weights.append(weight)
            i += 1
        return paths, weights
                
    def reachable_nodes_from(self, G):
        return {node: nx.descendants(G, node) | {node} for node in G.nodes()}

    def nodes_reaching(self, G):
        return {node: nx.ancestors(G, node) | {node} for node in G.nodes()}
    
    def min_cost_flow(self, G: nx.DiGraph, s, t, demands_attr='l', capacities_attr='u', costs_attr='c') -> tuple:
        bigNumber = 1 << 32

        flowNetwork = nx.DiGraph()
        flowNetwork.add_node(s, demand=-bigNumber)
        flowNetwork.add_node(t, demand=bigNumber)

        for v in G.nodes():
            if v != s and v != t:
                flowNetwork.add_node(v, demand=0)

        flowNetwork.add_edge(s, t, weight=0)

        counter = count(1)
        edgeMap = {}
        uid = "z" + str(id(G))

        for x, y in G.edges():
            z1 = uid + str(next(counter))
            z2 = uid + str(next(counter))
            edgeMap[(x, y)] = z1
            l = G[x][y][demands_attr]
            u = G[x][y][capacities_attr]
            c = G[x][y][costs_attr]
            flowNetwork.add_node(z1, demand=l)
            flowNetwork.add_node(z2, demand=-l)
            flowNetwork.add_edge(x, z1, weight=c, capacity=u)
            flowNetwork.add_edge(z1, z2, weight=0, capacity=u)
            flowNetwork.add_edge(z2, y, weight=0, capacity=u)

        try:
            flowCost, flowDictNet = nx.network_simplex(flowNetwork)
            flowDict = {node: {} for node in G.nodes()}
            for x, y in G.edges():
                flowDict[x][y] = flowDictNet[x][edgeMap[(x, y)]]
            return flowCost, flowDict
        except Exception:
            return None, None

    def compute_max_edge_antichain(self, G: nx.DiGraph, source, sink, get_antichain=False, weight_function=None) -> tuple:
        G_nx = nx.DiGraph()
        G_nx.add_nodes_from(G.nodes())
        demand = {}

        for u, v in G.edges():
            cost = 1 if u == source else 0
            edge_demand = int(u != source and v != sink)
            if weight_function:
                edge_demand = weight_function.get((u, v), 0)
            demand[(u, v)] = edge_demand
            G_nx.add_edge(u, v, l=edge_demand, u=(1 << 32), c=cost)

        minFlowCost, minFlow = self.min_cost_flow(G_nx, source, sink)

        def dfs_reachable(start, visited):
            stack = [start]
            while stack:
                u = stack.pop()
                if visited[u] != 0:
                    continue
                assert u != sink
                visited[u] = 1
                for v in G.successors(u):
                    if minFlow[u][v] > demand[(u, v)] and visited[v] == 0:
                        stack.append(v)
                for v in G.predecessors(u):
                    if visited[v] == 0:
                        stack.append(v)

        def dfs_saturating(start, visited, antichain):
            stack = [start]
            while stack:
                u = stack.pop()
                if visited[u] != 1:
                    continue
                visited[u] = 2
                for v in G.successors(u):
                    if minFlow[u][v] > demand[(u, v)]:
                        if visited[v] == 1:
                            stack.append(v)
                    elif minFlow[u][v] == demand[(u, v)] and demand[(u, v)] >= 1 and visited[v] == 0:
                        antichain.append((u, v))
                for v in G.predecessors(u):
                    if visited[v] == 1:
                        stack.append(v)

        if get_antichain:
            antichain = []
            visited = {node: 0 for node in G.nodes()}
            dfs_reachable(source, visited)
            dfs_saturating(source, visited, antichain)
            return minFlowCost, antichain

        return minFlowCost

    def _get_paths_to_fix_from_safe_lists(self, G, source, sink, safe_sequences) -> list:
        if not safe_sequences:
            return [], 0

        longest_safe_list = {}
        for i, safe_list in enumerate(safe_sequences):
            for edge in safe_list:
                if edge not in longest_safe_list or len(safe_sequences[longest_safe_list[edge]]) < len(safe_list):
                    longest_safe_list[edge] = i

        len_of_longest_safe_list = {
            edge: len(safe_sequences[longest_safe_list[edge]])
            for edge in longest_safe_list
        }

        _, edge_antichain = self.compute_max_edge_antichain(
            G, source, sink, get_antichain=True, weight_function=len_of_longest_safe_list
        )
        l = [safe_sequences[longest_safe_list[edge]] for edge in edge_antichain]
        
        length = 0
        for s in l:
            length += len(s)
            
        return l, length
    
    def solve(self, output: bool = False):
        """
        if self.G.number_of_edges() != 0 and not utils.check_st_graph(self.G):
            print("-------------------------------------")
            for u, v, data in self.G.edges(data=True):
                print(f"Edge: {u}->{v}, Data: {data}")
            print("-------------------------------------")
        """
            
        if self.G.number_of_edges() == 0:
            return None, None, None, None
        
        if self.greedy_width_solve:
            paths, weights = self.greedily_decompose_flow()
            return paths, None, None, weights

        aug_G = AbstractSourceSinkGraph(self.G)
        
        paths_to_fix, num_edges_fixed_with_one = self._get_paths_to_fix_from_safe_lists(aug_G, aug_G.source, aug_G.sink, self.safe_sequences)
        
        reachable = self.reachable_nodes_from(self.G)
        reaching  = self.nodes_reaching(self.G)
        k_start = max(1, len(paths_to_fix))
        num_edges_fixed_with_zero = 0
        for k in range(k_start, self.maximum_k + 1):
            myDecomp = kCFD.KCommonFlowDecomp(self.G, self.num_flows, k, self.flow_attr, self.subpath_constr, console=self.console)
            myDecomp.build_model()
            myDecomp.add_safe_sequences(paths_to_fix)
            num_edges_fixed_with_zero = myDecomp.fix_zero_edges(paths_to_fix, reachable, reaching)
            t0 = time.perf_counter()
            if myDecomp.solve_model():
                solve_time = time.perf_counter() - t0
                paths = myDecomp.get_model_paths()
                weights = myDecomp.get_model_weights()
                if self.console and output:
                    print(f"Found a solution with {k} distinct paths:\n" + myDecomp.get_model_solution())
                return paths, solve_time, num_edges_fixed_with_one + num_edges_fixed_with_zero, weights
        return None, None, None, None