from pyomo.environ import *

model = AbstractModel()

"""
Explanation: Max Flow
Summary
The goal of the maximum flow problem is to find the maximum flow possible in a network from some given source node to a given sink node.  Applications of the max flow problem include finding the maximum flow of orders through a job shop, the maximum flow of water through a storm sewer system, and the maximum flow of product through a product distribution system, among others. Schrijver (2002) note that the maximum flow problem was first formulated in 1954 by T. E. Harris and F. S. Ross as a simplified model of Soviet railway traffic flow.
A network is a directed graph, and the arc capacities, or upper bounds, are the only relevant parameters.  A network graph does not have to be symmetric: if an arc (v,w) is in the graph, the reverse arc (w,v) does not have to be in the graph. Further, parallel arcs are not allowed, and self-loops are not allowed. The source and the sink are distinct nodes in the network, but the sink may be unreachable from the source.
Problem Statement
The max flow problem can be formulated mathematically as a linear programming problem using the following model.  
Sets
$N$ = nodes in the network
 $A$ = network arcs
Parameters
$s$ = source node
 $t$ = sink node
 $c_{ij}$ = arc flow capacity, $\forall (i,j) \in A$
Variables
$f_{i,j}$ = arc flow, $\forall (i,j) \in A$
Objective
Maximize the flow into the sink nodes
 $\max \sum_{{i \mid (i,t) \in A}} c_{i,t} f_{i,t}$
Constraints
Enforce an upper limit on the flow across each arc
 $f_{i,j} \leq c_{i,j}$, $\forall (i,j) \in A$
Enforce flow through each node
 $\sum_{{i \mid (i,j) \in A}} f_{i,j} = \sum_{{i \mid (j,i) \in A}} f_{j,i}$, $\forall j \in N - {s,t}$
Flow lower bound
 $f_{i,j} \geq 0$, $\forall (i,j) \in A$
Pyomo Formulation
We begin by importing the Pyomo package and creating a model object:
End of explanation
"""


# Nodes in the network
model.N = Set()
# Network arcs
model.A = Set(within=model.N*model.N)

"""
Explanation: The sets $N$ and $A$ are declared abstractly using the Set component:
End of explanation
"""


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

"""
Explanation: Similarly, the model parameters are defined abstractly using the Param component:
End of explanation
"""


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

"""
Explanation: The within option is used in these parameter declarations to define expected properties of the parameters.  This information is used to perform error checks on the data that is used to initialize the parameter components.
The Var component is used to define the decision variables:
End of explanation
"""


# 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)

"""
Explanation: The within option is used to restrict the domain of the decision variables to the non-negative reals. This eliminates the need for explicit bound constraints for variables.
The Objective component is used to define the cost objective.  This component uses a rule function to construct the objective expression:
End of explanation
"""


# 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)

"""
Explanation: Similarly, rule functions are used to define constraint expressions in the Constraint component:
End of explanation
"""


!cat maxflow.py

"""
Explanation: Putting these declarations all together gives the following model:
End of explanation
"""


!cat maxflow.dat

"""
Explanation: Model Data
Since this is an abstract Pyomo model, the set and parameter values need to be provided to initialize the model.  The following data command file provides a synthetic data set:
End of explanation
"""


!pyomo solve --solver=glpk maxflow.py maxflow.dat

"""
Explanation: Set data is defined with the set command, and parameter data is defined with the param command.
This data set considers the problem of maximizing flow in a zoo. A panda is about to give birth at the zoo! Officials anticipate that attendance will skyrocket to see the new, adorable baby panda. There's one particular residential area called "Home" that is full of panda loving families and there's a fear that the increased number of people visiting the zoo will overload the public transportation system. It will be especially bad in the evening since the zoo closes about the same time as rush hour, so everyone will be trying to find space on the already crowded buses and subways. As a city planner, you were given a map of routes from the zoo to Home, along with the estimated number of families that could go on each route. Additionally, it was estimated that 16 families from Home will visit each day, and it's your task to figure out if this will overload the public transportation system, and, if it does, how could the system be improved?
Solution
Pyomo includes a pyomo command that automates the construction and optimization of models.  The GLPK solver can be used in this simple example:
End of explanation
"""


!cat results.yml

"""
Explanation: By default, the optimization results are stored in the file results.yml:
End of explanation
"""


