Constraint Programming OR-Tools finding Matrix Combinations

Question

As a similar post of mine Find all Combinations of a Matrix I am trying to find matrix combinations with entries $>0$ meaning for a matrix

\begin{bmatrix} 0 & 1 & 3 \\ 5 & 2 & 1 \\ 0 & 0 & 10 \end{bmatrix}

The following combinations: $(5,1,3); (5,2,1); (5,1,10);...$

I have decided to use constraint-programming with OR-Tools for this case to find all the solutions.

Furthermore, I want to use added constraints so that at most X rows should get used for a combination, e.g., max 2 rows $(5,1,3)$ but not (5,1,10).

At the moment I have the following model:

1. $x_{i,j}$ is a binary decision variable with $i$ being a row and $j$ being a column
2. $y_i$ a binary variable stating whether row $i$ is getting used for the combinations
3. I need exactly one row per column:

$$\sum_ix_{i,j} = 1 \ \forall j \in J$$

4. Choose only a matrix entity if it is $>0$ with $p_{i,j}$ stating the matrix elements: $$x_{i,j} \leq p_{i,j} \ \forall i,j$$

5. Help constraint to see whether a row is used for a combination by setting $y_i$ to $0$ or $1$ $$\sum_jx_{i,j} \leq y_i \cdot |J| \ \forall i \in I$$

6. Limit the use of rows for the combinations by a number of $q$. $$\sum_iy_i = q$$

1) (Solved) As of now, I am getting weird results when I define the variable $y_i$. By that, I mean that I am getting the same result multiple times.

2) Also, how to find which row is used in the combinations? As of now, the constraint 5 is inadequate, since y can still get the value 1 although the whole row is not chosen.

y = {}
       for i in range(matrixRows): # rows
          y[i] = model.NewBoolVar('vendor_%d' % (i))

If I comment it out the results are fine.

Here you can find the complete code with constraints 5 and 6 commented out (definition of $y_i$ included):

    from __future__ import print_function
    from ortools.sat.python import cp_model

    class SolutionPrinter(cp_model.CpSolverSolutionCallback):

        # def __init__(self,x,y,matrixRows,matrixColumns,matrix):
        def __init__(self,x,matrixRows,matrixColumns,matrix):
            cp_model.CpSolverSolutionCallback.__init__(self)
            self._x = x
            # self._y = y
            self._matrixRows = matrixRows
            self._matrixColumns = matrixColumns
            self._matrix = matrix
            self._solution_count = 0

        def OnSolutionCallback(self):
            self._solution_count += 1
            for j in range(self._matrixColumns):
                for i in range(self._matrixRows):
                    if self.Value(self._x[(i,j)]):
                        print('x: %d for x(%d,%d) and matrix value:%d' %(self.Value(self._x[(i,j)]),i,j,matrix[i][j]))
            # for v in self._x:
            # print('%s = %i' % (v, self.Value(v)), end=' ')
            print()

        def SolutionCount(self):
            return self._solution_count

    def main(matrix):
       model = cp_model.CpModel()
       matrixRows = len(matrix)
       matrixColumns = len(matrix[0])
      # define variables
       x = {}
       for i in range(matrixRows): # rows
         for j in range(matrixColumns):
           x[(i,j)] = model.NewBoolVar('company_%d,%d' % (i,j))

       y = {}
       for i in range(matrixRows): # rows
          y[i] = model.NewBoolVar('vendor_%d' % (i))

       # choose exactly one vendor per attribute
       for j in range(matrixColumns):
        model.Add(sum(x[(i,j)] for i in range(matrixRows)) == 1)

       # to choose value x(i,j), the parameter of the matrix has to be >0
       for i in range(matrixRows):
           for j in range(matrixColumns):
               model.Add(x[i,j] <= matrix[i][j])

       # # help constraint to count whether a vendor is in the mix
       # for i in range(matrixRows):
       #         model.Add(sum(x[i, j] for j in range(matrixRows)) <= y[i] *matrixColumns)
       #
       # # set the max vendors
       # model.Add(sum(y[i] for i in range(matrixRows)) <= 2)

       #call the solver and display the results
       solver = cp_model.CpSolver()
       # solution_printer = SolutionPrinter(x,y,matrixRows,matrixColumns,matrix)
       # the following is an object of the class SolutionPrinter
       solution_printer = SolutionPrinter(x,matrixRows,matrixColumns,matrix)
       status = solver.SearchForAllSolutions(model, solution_printer)
       print()
       print('Solutions found: %i' % solution_printer.SolutionCount())

    if __name__ == '__main__':
      matrix = [[0,1,3],[5,0,1],[0,0,10]]
      main(matrix)

Answer 1

You get the same result $2^{\operatorname{len}(y)}$ times because $y_i$ is not constrained.

You can see that by adding this to your callback:

print("y:", [self.Value(y) for y in self._y.values()])

x: 1 for x(1,0) and matrix value:5
x: 1 for x(0,1) and matrix value:1
x: 1 for x(1,2) and matrix value:1
y: [0, 0, 0]

x: 1 for x(1,0) and matrix value:5
x: 1 for x(0,1) and matrix value:1
x: 1 for x(1,2) and matrix value:1
y: [1, 0, 0]

...
y: [1, 0, 1]

...
y: [0, 0, 1]

...
y: [0, 1, 1]

...
y: [0, 1, 0]

...
y: [1, 1, 0]

...
y: [1, 1, 1]

PS: I would check matrix[i][j] != 0 before creating the boolean.

Edit:

model.Add(sum(x[i, j] for j in range(matrixRows)) > 0).OnlyEnforceIf(y[i])
model.Add(sum(x[i, j] for j in range(matrixRows)) == 0).OnlyEnforceIf(y[i].Not())