Skip to content

The Time Constrained Multi-Commodity Flow Problem in Liner Shipping

This examples illustrates an implementation of the problem presented in

Christian Vad Karsten, David Pisinger, Stefan Ropke, and Berit Dangaard Brouer, "The time constrained multi-commodity network flow problem and its application to liner shipping network design", Transportation Research Part E: Logistics and Transportation Review Volume 76, April 2015, Pages 122-138, https://doi.org/10.1016/j.tre.2015.01.005

The problem considered is the cargo-flow problem in a liner shipping network where there is a maximum limit that goods may be in transit. Goods are divided in origin-destination pairs with a given demand. The network consists of port calls with connecting edges given by rotations of the container ships, i.e, each leg sailed on a rotation corresponds to an edge in the network with the edge capacity equal to the ships capacity. Travel time on edges is given by the distance divided by the average sailing speed on the given rotation.

Port calls in the same port from different rotations are connected with transshipment edges to allow transshipment of goods from one ship to another. This is modeled as a complete subgraph between nodes in the same port. Load/unload edges connects origin and destination nodes with port call nodes and a direct forfeit edge between an origin-destination pair represents the option to leave the cargo on the dock.

The objective is to minimize the transshipment and the forfeiting cost while imposing a maximum transport time of goods. The problem is a variant of the multi-commodity flow problem with length constraints on the paths on which to flow.

Consider for each commodity kKk \in K the graph Gk(V,E)G^k(V,E) with vertices V=Vp{ok,dk}V = V_p \cup \{ o^k, d^k \} and edges E=EvEtEfElE = E_v \cup E_t \cup E_f \cup E_l where VpV_p is the set of port calls and (ok,dk)(o^k, d^k) is the origin and destination pair, EvE_v are the voyage edges between port calls, EtE_t are the transshipment edges, EfE_f are the forfeit edges, and ElE_l are the load/unload edges. Each edge ee have an associated cost cec_e and transshipment time tet_e. For voyage edges the cost is 0 and transshipment cost depend on vessel speed and distance. For transshipment edges the cost reflect the crane pricing and the transshipment time is the time for re-stowing the containers. The same is true for load/unload edges when moving containers from the dock to a ship. For forfeit edges a penalty is incurred for not delivering and the transshipment time is 0. Each commodity kKk \in K has a transshipment time limit TkT^k and a demand DkD^k associated.

The path based model is

minkKpPkλpks.t. pPkλpk=1kKkKDkαeλpkueeEv0λpk1kK,pPk \begin{aligned} \min{ } & \sum_{k \in K} \sum_{p \in P^k} \lambda_p^k \\ \text{s.t. } & \sum_{p \in P^k} \lambda_p^k = 1 & & k \in K \\ & \sum_{k \in K} D^k \alpha_e \lambda_p^k \leq u_e & & e \in E_v \\ & 0 \leq \lambda_p^k \leq 1 && k \in K, p \in P^k \end{aligned}

where pPkp \in P^k are time feasible paths with source oko^k and sink dkd^k in graphs Gk(V,E)G^k(V,E) with accumulated edge cost cpk=eEce(sp:s=e1)c_p^k = \sum_{e \in E} c_e\Big ( \sum_{ s \in p: s = e} \mathbf{1} \Big ) and αe=sp:s=e1\alpha_e = \sum_{ s \in p: s = e} \mathbf{1} indicates usage of edge ee on pp.

The equivalent network model is

minkKpPkeEcekDkxeks.t. eδ+(ok)xek=1kKkKDkxekueeEvxek=pPk(sp:s=e1)λpkkK,eE1pPkλpk1kK0λpk1kK,pPk0xek1kK,eE \begin{aligned} \min{ } & \sum_{k \in K} \sum_{p \in P^k} \sum_{e \in E} c_e^k D^k x_e^k \\ \text{s.t. } & \sum_{e \in \delta^+(o^k)} x_e^k = 1 & & k \in K \\ & \sum_{k \in K} D^k x_e^k \leq u_e & & e \in E_v \\ & x_e^k = \sum_{p \in P^k} \Big ( \sum_{ s \in p: s = e} \mathbf{1} \Big ) \lambda_p^k & & k \in K, e \in E\\ & 1 \leq \sum_{p \in P^k} \lambda_p^k \leq 1 && k \in K \\ & 0 \leq \lambda_p^k \leq 1 && k \in K, p \in P^k \\ & 0 \leq x_e^k \leq 1 && k \in K, e \in E \end{aligned}

where pPkp \in P^k are paths with source oko^k and sink dkd^k in graphs Gk(V,E)G^k(V,E) subject to the time constraint modeled as a disposable resource with edge consumption tet_e and vertex bounds [0,Tk][0,T^k] for all iVi \in V.

The model to implement is

minkKpPkeEcekDkxeks.t. eδ+(ok)xek=1kKkKDkxekueeEv0xek1kK,eE \begin{aligned} \min{ } & \sum_{k \in K} \sum_{p \in P^k} \sum_{e \in E} c_e^k D^k x_e^k \\ \text{s.t. } & \sum_{e \in \delta^+(o^k)} x_e^k = 1 & & k \in K \\ & \sum_{k \in K} D^k x_e^k \leq u_e & & e \in E_v \\ & 0 \leq x_e^k \leq 1 && k \in K, e \in E \end{aligned}

with the graph and resource constraints described above. See the modeling section for details.

This examples depends on networkx. Install with

pip install networkx

The code is

import networkx
from flowty import Model, xsum, LinExpr
from flowty.datasets import linerlib

# data = linerlib.fetch_linerlib(instance="Mediterranean")
# network = linerlib.fetch_linerlib_rotations(instance="Med_base_best")

data = linerlib.fetch_linerlib(instance="Baltic")
network = linerlib.fetch_linerlib_rotations(instance="Baltic_best_base")

name, _, _, _, _ = data["instance"]

builder = linerlib.GraphBuilder(data, network)

# build the graph
g = networkx.DiGraph()

# port call, origin, destination nodes
g.add_nodes_from(
    [
        (n, {"index": g.number_of_nodes() + i})
        for i, n in enumerate(builder.portCallNodes())
    ]
)
g.add_nodes_from(
    [
        (n, {"index": g.number_of_nodes() + i})
        for i, n in enumerate(builder.originNodes())
    ]
)
g.add_nodes_from(
    [
        (n, {"index": g.number_of_nodes() + i})
        for i, n in enumerate(builder.destinationNodes())
    ]
)

# voyage, transit, load, forfeit edges
voyageEdges = builder.voyageEdges()
g.add_edges_from(
    [
        (n[0], n[1], {"index": g.number_of_edges() + i})
        for i, n in enumerate(voyageEdges)
    ]
)
g.add_edges_from(
    [
        (n[0], n[1], {"index": g.number_of_edges() + i})
        for i, n in enumerate(builder.transitEdges())
    ]
)
g.add_edges_from(
    [
        (n[0], n[1], {"index": g.number_of_edges() + i})
        for i, n in enumerate(builder.loadEdges())
    ]
)
g.add_edges_from(
    [
        (n[0], n[1], {"index": g.number_of_edges() + i})
        for i, n in enumerate(builder.forfeitEdges())
    ]
)

# model building
m = Model()

# number of subproblems
k = len(builder.demand["Destination"])

# for keeping track of voyage edge variables
voyageEdgeIds = [g.edges[edge[0], edge[1]] for edge in voyageEdges]
voyageEdgeVarsIds = [{} for j in range(len(voyageEdgeIds))]

# resource constrained graphs
gs = []
for i in range(k):
    origin = builder.demand["Origin"][i]
    dest = builder.demand["Destination"][i]
    vs = builder.portCallNodes() + [f"O{i}_{origin}", f"D{i}_{dest}"]
    es = g.edges(nbunch=vs)

    source = g.nodes[f"O{i}_{origin}"]["index"]
    sink = g.nodes[f"D{i}_{dest}"]["index"]
    obj = [
        builder.demand["FFEPerWeek"][i] * builder.cost[g.edges[e]["index"]] for e in es
    ]
    edges = [(g.nodes[e[0]]["index"], g.nodes[e[1]]["index"]) for e in es]
    gk = m.addGraph(
        obj=obj,
        edges=edges,
        source=source,
        sink=sink,
        L=1,
        U=1,
        type="C",
        names=f"x_{i}",
    )

    time = [builder.travelTime[g.edges[e]["index"]] for e in es]
    m.addResourceDisposable(
        graph=gk,
        consumptionType="E",
        weight=time,
        boundsType="V",
        lb=0,
        ub=builder.demand["TransitTime"][i],
        name=f"time_{i}",
    )
    gs.append(gk)

    # keep track of voyageEdges for constraints later
    tmp = list([eid["index"] for eid in voyageEdgeIds])
    for j, e in enumerate(es):
        for h, eid in enumerate(tmp):
            if g.edges[e]["index"] == eid:
                voyageEdgeVarsIds[eid][i] = j
                tmp.remove(eid)
                break

# graph vars
vars = [gs[i].vars for i in range(k)]

# sum_( i,j \in delta+(o^k)) x_ijk = 1 , forall k
for i in range(k):
    source = gs[i].source
    m.addConstr(xsum((1, x) for x in vars[i] if source == x.source) == 1)

for j, ks in enumerate(voyageEdgeVarsIds):
    expr = LinExpr()
    for i, h in ks.items():
        x = vars[i][h]
        expr.addTerm(builder.demand["FFEPerWeek"][i], x)
    m.addConstr(expr <= builder.capacity[j])

status = m.optimize()

print(f"ObjectiveValue {m.objectiveValue}")

# get the variables
for var in m.vars:
    if var.x > 0:
        print(f"{var.name} = {var.x}")