pythonalgorithmmathematical-optimizationoperations-research

Optimal combination of linked-buckets


Let's say I have the following (always binary) options:

import numpy as np
a=np.array([1, 1, 0, 0, 1, 1, 1])
b=np.array([1, 1, 0, 0, 1, 0, 1])
c=np.array([1, 0, 0, 1, 0, 0, 0])
d=np.array([1, 0, 1, 1, 0, 0, 0])

And I want to find the optimal combination of the above that get's me to at least, with minimal above:

req = np.array([50,50,20,20,100,40,10])

For example:

final = X1*a + X2*b + X3*c + X4*d
  1. Does this map to a known operational research problem? Or does it fall under mathematical programming?
  2. Is this NP-hard, or exactly solveable in a reasonable amount of time (I've assumed it's combinatorally impossible to solve exactly)
  3. Are there know solutions to this?

Note: The actual length of arrays are longer - think ~50, and the number of options are ~20 My current research has led me to some combination of the assignment problem and knapsack, but not too sure.


Solution

  • It's a covering problem, easily solvable using an integer program solver (I used OR-Tools below). If the X variables can be fractional, substitute NumVar for IntVar. If the X variables are 0--1, substitute BoolVar.

    import numpy as np
    
    a = np.array([1, 1, 0, 0, 1, 1, 1])
    b = np.array([1, 1, 0, 0, 1, 0, 1])
    c = np.array([1, 0, 0, 1, 0, 0, 0])
    d = np.array([1, 0, 1, 1, 0, 0, 0])
    opt = [a, b, c, d]
    req = np.array([50, 50, 20, 20, 100, 40, 10])
    
    
    from ortools.linear_solver import pywraplp
    
    solver = pywraplp.Solver.CreateSolver("SCIP")
    x = [solver.IntVar(0, solver.infinity(), "x{}".format(i)) for i in range(len(opt))]
    extra = [solver.NumVar(0, solver.infinity(), "y{}".format(j)) for j in range(len(req))]
    for j, (req_j, extra_j) in enumerate(zip(req, extra)):
        solver.Add(extra_j == sum(opt_i[j] * x_i for (opt_i, x_i) in zip(opt, x)) - req_j)
    solver.Minimize(sum(extra))
    status = solver.Solve()
    if status == pywraplp.Solver.OPTIMAL:
        print("Solution:")
        print("Objective value =", solver.Objective().Value())
        for i, x_i in enumerate(x):
            print("x{} = {}".format(i, x[i].solution_value()))
    else:
        print("The problem does not have an optimal solution.")
    

    Output:

    Solution:
    Objective value = 210.0
    x0 = 40.0
    x1 = 60.0
    x2 = -0.0
    x3 = 20.0