Using networkx predecessors in Pyomo initialize method

Question

I am working on the directed graph by using the Networkx package and what I need is to use its predecessors' method on an optimization model. Let's say, there exists a directed graph with just $12$ nodes which is defined as follows:

G = nx.DiGraph(
    [(1, 2), (1, 3), 
     (2, 4), (2, 5), (2, 6), 
     (3, 6), (3, 7), 
     (4, 9), 
     (5, 9), 
     (6, 8),
     (7, 11),
     (8, 9), (8, 10),
     (9, 12),
     (10, 12),
     (11, 12),
     ])

Now, I would like to use its predecessors in the constraints:

$$ \sum_{m} m.x_{i,m} \leq \sum_{m} m.x_{j,m} \quad \forall (i,j) \in \text{Predecessors}$$

What I try is something like:

predecessors = numpy.array(G.predecessors)
m.pred= Set(initialize = m.I * m.I, dimen=2, filter=lambda m, i, j: (i,j) in predecessors)

By that, the model does not produce any throw, but the results are incorrect. I guess it came from the definition of the above set (!) and am wondering if, someone can guide me to fix that.

As an attempt, I tried to run the problem in DoCplex, but the issue on the $3$th constraint is the same. I think(?) this should be more of a Python issue than a solver issue.

number_of_nodes = 12
S = [[1, 2], [1, 3], 
     [2, 4], [2, 5], [2, 6], 
     [3, 6], [3, 7], 
     [4, 9], 
     [5, 9], 
     [6, 8],
     [7, 11],
     [8, 9], [8, 10],
     [9, 12],
     [10, 12],
     [11, 10],
]

I = range(number_of_nodes)
B = range(number_of_station)

x = {}
for i in I:
  for b in B:
    x[(i,b)] = mdl.binary_var()

the first form:
for i in S:
  for j in S:
    if (i,j) in enumerate(S):
      mdl.add_constraint( mdl.sum(x[(i,b)]*b for b in B) <= mdl.sum(x[(j,b)]*b for b in B ) ) 

the second form:
for i in S:
  for j in i:
    if (i,j) in enumerate(S):
      mdl.add_constraint( mdl.sum(x[(i,b)]*b for b in B) <= mdl.sum(x[(j,b)]*b for b in B ) )

Unfortunately, in both form the solution is incorrect.

The correct solution:
Objective = 20

----------------------
                 1           2           3           4

    1        1.000
    2        1.000
    3                    1.000
    4        1.000
    5        1.000
    6                                1.000
    7                    1.000
    8                                1.000
    9                                            1.000
    10                                           1.000
    11                   1.000
    12                                           1.000
  

Answer 1

As you already guessed, these are rather coding issues than solver issues.

Regarding the pyomo model: Both m.pred and m.PAIRS in your linked MWE are empty. DiGraph's predecessor method returns an iterator and expects the node n as argument whose predecessors you'd like to iterate over. That's why

predecessors = numpy.array(G.predecessors)

is just a numpy array of a method and not an iterator. So it doesn't make sense in your case. Since you can initialize Pyomo Sets from any Iterable, you can use a simple list instead of an iterator:

predecessors = list(list(G.predecessors(node)) for node in G.nodes())

Then, predecessors[i-1] gives you all the direct predecessor nodes of the node $i$ as Python uses zero-based indexing and your graph nodes start with 1. Alternatively, you can easily create a dict such that predecessors[i, j] gives you the direct predecessor nodes of the nodes $i$ and $j$, i.e. the edge $(i, j)$:

predecessors = {(i,j): list(G.predecessors(i)) + list(G.predecessors(j)) for (i,j) in S}

Regarding the DocPlex model it's not really clear what if (i, j) in enumerate(S) is expected to do. enumerate is just an iterator (not a list!) that gives you the index of each edge in S and the edge stored in a python list.

Depending on what "$\forall(i, j) \in \text{Precedessors}$" means (it's a bit vague, IMO), the above tips may be helpful in getting you back on track.