Min Cost Flow Slotting › Implementation in OR-Tools

from ortools.graph.python import min_cost_flow


def solve_minimum_cost_flow(
    start_nodes, end_nodes, capacities, unit_costs, supplies
):
    """
    Solve a minimum cost flow problem using OR-Tools.

    This function translates our mathematical formulation into OR-Tools API calls:
    - Edges (i,j) with capacity u_ij and cost c_ij → add_arc_with_capacity_and_unit_cost()
    - Node supplies b_i → set_node_supply()
    - The solver automatically handles flow conservation, capacity constraints, and optimization

    Parameters:
    -----------
    start_nodes : list
        List of starting nodes for each arc (the "i" in edge (i,j))
    end_nodes : list
        List of ending nodes for each arc (the "j" in edge (i,j))
    capacities : list
        List of capacities u_ij for each arc (upper bounds on flow)
    unit_costs : list
        List of unit costs c_ij for each arc (cost per unit of flow)
    supplies : list
        List of supplies b_i for each node (positive for sources, negative for sinks, zero for transshipment)

    Returns:
    --------
    dict : Dictionary with solution information including optimal flows and total cost
    """
    # Create the solver instance
    # This initializes the internal data structures that will represent our network
    solver = min_cost_flow.SimpleMinCostFlow()

    # Step 1: Add all edges to the network
    # This corresponds to defining the set E in our graph G = (V, E)
    # For each edge, we specify: start node, end node, capacity (u_ij), and unit cost (c_ij)
    for i in range(len(start_nodes)):
        solver.add_arc_with_capacity_and_unit_cost(
            start_nodes[i], end_nodes[i], capacities[i], unit_costs[i]
        )

    # Step 2: Set node supplies/demands
    # This corresponds to setting the values b_i for each node i ∈ V
    # The solver will use these to enforce flow conservation constraints
    for i in range(len(supplies)):
        solver.set_node_supply(i, supplies[i])

    # Step 3: Solve the optimization problem
    # The solver automatically:
    # - Enforces flow conservation: Σx_out - Σx_in = b_i for each node
    # - Enforces capacity constraints: 0 ≤ x_ij ≤ u_ij for each edge
    # - Minimizes the objective: min Σ c_ij * x_ij
    status = solver.solve()

    # Step 4: Extract and return the solution
    if status == solver.OPTIMAL:
        # The solver found an optimal solution
        solution = {
            "status": "OPTIMAL",
            "total_cost": solver.optimal_cost(),  # The value of our minimized objective function
            "flows": {},  # Dictionary mapping (start, end) → flow amount x_ij
            "arc_info": [],  # List of detailed information for each arc
        }

        # Extract the optimal flow value for each edge
        # These are the values of our decision variables x_ij that minimize cost
        for i in range(solver.num_arcs()):
            tail = solver.tail(i)  # Start node of this arc
            head = solver.head(i)  # End node of this arc
            flow = solver.flow(i)  # Optimal flow x_ij on this arc
            capacity = solver.capacity(i)  # Capacity u_ij of this arc
            unit_cost = solver.unit_cost(i)  # Unit cost c_ij of this arc
            arc_cost = flow * unit_cost  # Total cost contribution: c_ij * x_ij

            # Store the flow in our solution dictionary
            solution["flows"][(tail, head)] = flow
            solution["arc_info"].append(
                {
                    "arc": (tail, head),
                    "flow": flow,
                    "capacity": capacity,
                    "unit_cost": unit_cost,
                    "total_cost": arc_cost,
                }
            )

        return solution
    else:
        # The problem is infeasible (no solution exists) or the solver encountered an error
        return {
            "status": "INFEASIBLE"
            if status == solver.INFEASIBLE
            else "UNKNOWN",
            "total_cost": None,
            "flows": {},
            "arc_info": [],
        }


# Example: Supply chain optimization
def supply_chain_example():
    """
    Example of supply chain optimization using MCF.

    This implements the supply chain problem we worked through in the Example section:
    - Factory (node 0) produces 10 units (supply = 10)
    - Warehouse (node 1) and Distribution (node 2) are transshipment nodes (supply = 0)
    - Market (node 3) requires 10 units (supply = -10, which represents demand)

    The network has 5 edges representing different shipping routes, each with
    capacity limits and transportation costs. The solver will find the cheapest
    way to route all 10 units from factory to market.
    """

    # Define the network structure
    # Each list has the same length, with index i representing the i-th edge
    start_nodes = [0, 0, 1, 1, 2]  # From nodes: edges start at these nodes
    end_nodes = [1, 2, 2, 3, 3]  # To nodes: edges end at these nodes
    # Together, these define edges: (0,1), (0,2), (1,2), (1,3), (2,3)

    # Capacity constraints: u_ij for each edge
    # These are the upper bounds that prevent overloading any route
    capacities = [15, 8, 5, 10, 12]  # Maximum flow allowed on each edge

    # Unit costs: c_ij for each edge
    # These determine the cost per unit of flow, used in the objective function
    unit_costs = [4, 6, 1, 2, 3]  # Cost per unit of flow on each edge

    # Node supplies: b_i for each node
    # Positive = supply (source), negative = demand (sink), zero = transshipment
    supplies = [10, 0, 0, -10]  # Factory supplies 10, Market demands 10

    # Solve the optimization problem
    # This function call translates our problem into the mathematical formulation
    # and uses OR-Tools to find the optimal flow values x_ij
    solution = solve_minimum_cost_flow(
        start_nodes, end_nodes, capacities, unit_costs, supplies
    )

    return solution


# Run the example
# This solves the supply chain problem we analyzed mathematically earlier
supply_chain_solution = supply_chain_example()

def time_slotting_example():
    """
    Example of resource allocation across time slots using MCF.

    This demonstrates how MCF can model slotting problems where resources
    must be allocated across time periods. The network structure is:

    Source → Resources → Time Slots → Tasks → Sink

    This creates a multi-stage flow where:
    1. Resources are allocated from a source
    2. Resources are assigned to time slots (with costs that may vary by time)
    3. Time slots are connected to tasks (with costs that may vary by task priority)
    4. Tasks flow to a sink (representing completed work)

    The solver finds the optimal allocation that minimizes total cost while
    respecting capacity constraints at each stage.
    """

    # Define the problem dimensions
    time_slots = 3  # We have 3 time periods to allocate resources
    resources = 2  # We have 2 resources (e.g., workers, machines)
    tasks = 4  # We need to complete 4 tasks

    # Build the network structure layer by layer
    # This creates a multi-stage network: Source → Resources → Time Slots → Tasks → Sink
    start_nodes = []
    end_nodes = []
    capacities = []
    unit_costs = []

    # Stage 1: Source to Resources
    # The source node (node 0) represents the total available capacity
    # We connect it to each resource, allowing flow to be allocated to different resources
    source = 0
    for r in range(resources):
        start_nodes.append(source)
        end_nodes.append(r + 1)  # Resource nodes are numbered 1, 2, ...
        capacities.append(
            10
        )  # Each resource has capacity to handle 10 units total
        unit_costs.append(
            0
        )  # No cost for allocating resources (costs occur when using them)

    # Stage 2: Resources to Time Slots
    # Each resource can be assigned to work in different time slots
    # This creates edges from each resource to each time slot
    # Node numbering: resources are 1..resources, time slots are (resources+1)..(resources+time_slots)
    for r in range(resources):
        for t in range(time_slots):
            start_nodes.append(r + 1)  # From resource r
            end_nodes.append(resources + 1 + t)  # To time slot t
            capacities.append(
                5
            )  # Each resource can work 5 units per time slot (capacity constraint)
            unit_costs.append(
                1
            )  # Base cost for a resource working in a time slot

    # Stage 3: Time Slots to Tasks
    # Each time slot can contribute to completing different tasks
    # This creates edges from each time slot to each task
    # Node numbering: tasks are (resources+time_slots+1)..(resources+time_slots+tasks)
    for t in range(time_slots):
        for task in range(tasks):
            start_nodes.append(resources + 1 + t)  # From time slot t
            end_nodes.append(resources + time_slots + 1 + task)  # To task
            capacities.append(3)  # Each time slot can handle 3 units per task
            unit_costs.append(
                t + 1
            )  # Cost increases with later time slots (urgency penalty)
            # This models that completing work earlier is preferred

    # Stage 4: Tasks to Sink
    # The sink node represents completed work
    # All tasks must flow to the sink, representing that all tasks must be completed
    sink = resources + time_slots + tasks + 1
    for task in range(tasks):
        start_nodes.append(resources + time_slots + 1 + task)  # From task
        end_nodes.append(sink)  # To sink
        capacities.append(8)  # Each task can accept up to 8 units of work
        unit_costs.append(
            0
        )  # No cost for completing tasks (costs are in the process)

    # Define node supplies/demands
    # This enforces that:
    # - Source supplies total capacity (20 units)
    # - All intermediate nodes are transshipment (supply = 0)
    # - Sink demands the same amount (20 units) to ensure all work is completed
    supplies = [20] + [0] * (resources + time_slots + tasks) + [-20]

    # Solve the optimization problem
    # The solver will find the optimal way to route flow through all stages
    # while minimizing total cost and respecting all capacity constraints
    solution = solve_minimum_cost_flow(
        start_nodes, end_nodes, capacities, unit_costs, supplies
    )

    # Analyze the solution by aggregating flows at each stage
    # This helps us understand how the solver allocated resources
    resource_usage = {}  # How much each resource is used
    time_slot_usage = {}  # How much work is done in each time slot
    task_usage = {}  # How much work is allocated to each task

    # Parse the solution to extract flow information
    # We examine each edge in the solution and aggregate flows by stage
    if solution["status"] == "OPTIMAL":
        for arc_info in solution["arc_info"]:
            arc = arc_info["arc"]
            flow = arc_info["flow"]  # The optimal flow x_ij on this edge

            if flow > 0:  # Only consider edges that are actually used
                start, end = arc

                # Stage 1: Source to Resources
                # Track how much flow is allocated to each resource
                if start == 0:  # Edge starts at source
                    resource_usage[end] = resource_usage.get(end, 0) + flow

                # Stage 2: Resources to Time Slots
                # Track how much work each resource does in each time slot
                elif 1 <= start <= resources:  # Edge starts at a resource
                    time_slot = end - resources - 1  # Extract time slot index
                    time_slot_usage[time_slot] = (
                        time_slot_usage.get(time_slot, 0) + flow
                    )

                # Stage 3: Time Slots to Tasks
                # Track how much work each time slot contributes to each task
                elif (
                    resources + 1 <= start <= resources + time_slots
                ):  # Edge starts at a time slot
                    task = (
                        end - resources - time_slots - 1
                    )  # Extract task index
                    task_usage[task] = task_usage.get(task, 0) + flow

    return (
        solution,
        resource_usage,
        time_slot_usage,
        task_usage,
        resources,
        time_slots,
        tasks,
    )


# Run the time slotting example
(
    time_slotting_solution,
    resource_usage,
    time_slot_usage,
    task_usage,
    num_resources,
    num_time_slots,
    num_tasks,
) = time_slotting_example()

def manual_network_construction():
    """
    Manually construct and solve a simple MCF problem using direct solver API.

    This demonstrates the direct mapping between mathematical concepts and code:
    - Each edge (i,j) → add_arc_with_capacity_and_unit_cost(i, j, u_ij, c_ij)
    - Each node supply b_i → set_node_supply(i, b_i)
    - Query results → flow(i), optimal_cost()

    This approach is useful when you want fine-grained control over the problem
    construction or when building networks incrementally.
    """

    # Initialize the solver
    # This creates an empty network that we'll build edge by edge
    solver = min_cost_flow.SimpleMinCostFlow()

    # Step 1: Add edges (arcs) to the network
    # Each call adds one edge (i,j) with capacity u_ij and unit cost c_ij
    # This directly corresponds to defining the set E in our graph G = (V, E)

    solver.add_arc_with_capacity_and_unit_cost(
        0, 1, 10, 2
    )  # Edge (0,1): u=10, c=2
    solver.add_arc_with_capacity_and_unit_cost(
        0, 2, 8, 4
    )  # Edge (0,2): u=8, c=4
    solver.add_arc_with_capacity_and_unit_cost(
        1, 2, 5, 1
    )  # Edge (1,2): u=5, c=1
    solver.add_arc_with_capacity_and_unit_cost(
        1, 3, 12, 3
    )  # Edge (1,3): u=12, c=3
    solver.add_arc_with_capacity_and_unit_cost(
        2, 3, 15, 2
    )  # Edge (2,3): u=15, c=2

    # Step 2: Set node supplies/demands
    # This sets the values b_i for each node, which the solver uses to enforce
    # flow conservation: Σx_out - Σx_in = b_i for each node i

    solver.set_node_supply(0, 15)  # Node 0: supply node (b_0 = 15)
    solver.set_node_supply(1, 0)  # Node 1: transshipment (b_1 = 0)
    solver.set_node_supply(2, 0)  # Node 2: transshipment (b_2 = 0)
    solver.set_node_supply(3, -15)  # Node 3: demand node (b_3 = -15)

    # Step 3: Solve the optimization problem
    # The solver automatically:
    # - Enforces flow conservation at each node
    # - Enforces capacity constraints on each edge
    # - Minimizes total cost: min Σ c_ij * x_ij
    status = solver.solve()

    # Extract solution information
    solution_data = {
        "status": status,
        "optimal": status == solver.OPTIMAL,
        "cost": solver.optimal_cost() if status == solver.OPTIMAL else None,
        "arcs": [],
    }

    if status == solver.OPTIMAL:
        # Extract the solution: the optimal values of our decision variables x_ij
        # The solver stores arcs in the order they were added
        for i in range(solver.num_arcs()):
            tail = solver.tail(i)  # Start node of arc i
            head = solver.head(i)  # End node of arc i
            flow = solver.flow(i)  # Optimal flow x_ij on this arc
            capacity = solver.capacity(i)  # Capacity u_ij
            unit_cost = solver.unit_cost(i)  # Unit cost c_ij
            total_cost = flow * unit_cost  # Cost contribution: c_ij * x_ij

            solution_data["arcs"].append(
                {
                    "tail": tail,
                    "head": head,
                    "flow": flow,
                    "capacity": capacity,
                    "unit_cost": unit_cost,
                    "total_cost": total_cost,
                }
            )

    return solution_data


# Run manual example
# This demonstrates the direct API usage, which is useful for understanding
# how the mathematical formulation maps to code
manual_solution = manual_network_construction()