Improve Network Flow computation with additional Pyomo constraints based on previous computation

Question

I am performing several max-flow computations on extremely similar networks and I seek to improve their execution time based on information from previous computations (As there is a lot of repetitive computation). The code below is from the following repository + some additional details that I've added: https://github.com/Pyomo/PyomoGallery/blob/master/maxflow/maxflow.py

I am interested in the solutions for computations which are equal to a predefined_sum parameter. In addition I've added the additional_rule which aims to apply a constraint that the maximum flow in the sink node should be equal to predefined_sum

However, this seems to slow down the computation and does not improve the running time. Does someone have information about why this happens? Additionally, please let me know if you know some additional way to speed up the computation. I was also advised that the CPLEX API offers more functionality for such improvements so if someone knows some specific resources that can be used would be very helpful.

from pyomo.environ import *

model = AbstractModel()
#Specified flow
model.predefined_sum = pyo.Param(within=NonNegativeReals)
# Nodes in the network
model.N = Set()
# Network arcs
model.A = Set(within=model.N*model.N)

# Source node
model.s = Param(within=model.N)
# Sink node
model.t = Param(within=model.N)
# Flow capacity limits
model.c = Param(model.A)

# The flow over each arc
model.f = Var(model.A, within=NonNegativeReals)

# Maximize the flow into the sink nodes
def total_rule(model):
    return sum(model.f[i,j] for (i, j) in model.A if j == value(model.t))
model.total = Objective(rule=total_rule, sense=maximize)

# Enforce an upper limit on the flow across each arc
def limit_rule(model, i, j):
    return model.f[i,j] <= model.c[i, j]
model.limit = Constraint(model.A, rule=limit_rule)

# Enforce flow through each node
def flow_rule(model, k):
    if k == value(model.s) or k == value(model.t):
        return Constraint.Skip
    inFlow  = sum(model.f[i,j] for (i,j) in model.A if j == k)
    outFlow = sum(model.f[i,j] for (i,j) in model.A if i == k)
    return inFlow == outFlow
model.flow = Constraint(model.N, rule=flow_rule)

def additional_rule(model):
    return sum(model.f[i,j] for (i, j) in model.A if j == value(model.t)) == model.predefined_sum
model.additional_rule = pyo.Constraint(rule=additional_rule)

```

Answer 1

Honestly, I do not know how to embed previous computations results in a new execution when considering a mathematical model. However, you can reuse previous computations in a new iteration, when considering an "analytical approach", i.e. one not based on the use of solvers. The strategy relies, in principle, on any algorithm that makes use of $s$-$t$-paths in $G$, where $s$ and $t$ are the source and sink nodes, respectively. For didactical purposes, we will take the Ford-Fulkerson algorithm, where we will try to reuse the previous computations whenever calculating a new $s$-$t$-path in $G$. For a better explanation, let's take the Ford-Fulkerson algorithm given below.

• Inputs: Network $G=(V,E)$ with flow capacity $c$, a source node $s$, and a sink node $t$;
• Output: Compute a flow $f$ from $s$ to $t$ of maximum value.

1. $f(u,v)\leftarrow 0$ for all edges $(u,v)$;
2. Let $G_{f}(V,E_{f})$ be the residual network with capacity $c_{f}(u,v) = c(u,v) - f(u,v)$
3. While there is a path $p$ from $s$ to $t$ in $G_{f}$, such that $c_{f}(u,v)>0$ for all edges $(u,v) \in p$:
   1. Find $c_{f}(p)=\min\{c_{f}(u,v):(u,v)\in p\}$;
   2. For each edge $(u,v)\in p$:
      1. $f(u,v)\leftarrow f(u,v)+c_{f}(p)$ (Send flow along the path);
      2. $f(v,u)\leftarrow f(v,u)-c_{f}(p)$ (The flow might be "returned" later);

Let's say that the algorithm above is executed $N$ times, and each Max-Flow execution $i$ caches the generated $s$-$t$-paths in a pool of paths, named as, $P_i$, i.e., the $i$-th execution saves the paths generated at step 3 in the pool $P_i$.

In a given Max-Flow execution $j$, we may reduce the number of iterations made at step 3, by making look-ups in the pools $\{P_i: i < j\}$, for any path in which all the nodes belong to the $j$-th execution input graph, that is, $\{ p \in P_i: V(p) \subseteq V_j \wedge i < j\}$. Of course, a smart strategy for reducing the look-up effort may be required, maybe a trie tree could suit this case.

Not that, the described reusing strategy is not limited to the Ford-Fulkerson algorithm, but to algorithms that rely on the use of $s$-$t$-paths.

Case any point is not clear, please, let me know. Regards.

Answer 2

Try using pyomo quicksum instead of sum for defining the additional_rule. Quicksum improves efficiency with scale. Also please check if you are using latest and compatible Python & Pyomo version.
Also check if in the additional_rule using $\le$ rather than $==$ makes any difference, since objective is being maximized. Tight constraints sometimes take more time.