Trouble with mixed integer linear programming

asked 2018-03-10 01:53:02 +0100

JRHales
41 ●4 ●5 ●10

updated 2018-03-10 16:25:06 +0100

I'm having trouble with sage's mixed linear programming tool. I entered the following code trying to maximize a variable with the given constraint that multiplication by a matrix would result in the zero vector.

picture of inputs is here

sage: p = MixedIntegerLinearProgram(maximization = True, solver = "GLPK")
sage: x = p.new_variable(integer = True, nonnegative = True)
sage: p.add_constraint(x*spoof_matrix == 0)
sage: p.set_objective(x[25])
sage: p.solve()

Where spoof_matrix is a 56x34 matrix that has the specific constraints I'm working with.

When I use the p.solve() it gives the error in the picture.

I know that there is at least one non-trivial solution since the test vector seen in the picture solves the constraints.

I don't know why sage isn't solving the optimization problem, and any help would be really great!

Also, if it's possible to not only get the maximum value but also obtain the vector that gives the maximum value, that would be really awesome. Thanks!!

Edit: People have asked for the code for how spoof matrix is created, using these functions the spoof matrix would be given by the following line.

create_spoof_matrix(reduce_n_times(25,generate_possible_factors_list(3),primes_first_n(168),168),168)

def spoof(p,a):
    if p*a == 0:
        return 1
    if p!=1 and a != 1:
        value = ((1/p)^(a+1)-1)/(1/p-1)
    elif p ==1:
        value = (a+1)
    elif a ==1:
        value = (1+1/p)
    return value

def reduce_possible_primes(possible_factors,number_of_primes):
    new_prime_list = []
    for x in possible_factors:
        for p in primes_first_n(number_of_primes):
            if x % p == 0:
                if p not in new_prime_list:
                    new_prime_list.append(p)
    for x in xrange(0,len(new_prime_list)): #this is just ordering the primes from least to greatest.
        for y in xrange(0,len(new_prime_list)):
            if new_prime_list[x] < new_prime_list[y]:
                cache = new_prime_list[x]
                new_prime_list[x] = new_prime_list[y]
                new_prime_list[y] = cache
    return new_prime_list


def generate_possible_factors_list(size_of_factors):
    possible_factors = []
    for x in xrange(-1*(10^size_of_factors)+1,10^size_of_factors,2):
        if x != -1:
            numerator_factor = list(factor(numerator(spoof(x,2))))
            if numerator_factor[-1][0] < 10^size_of_factors:
                possible_factors.append(x)
    return possible_factors



def reduce_possible_factors(possible_factors,possible_primes):
    for x in possible_factors:
        numerator_factor_bank = list(factor(numerator(spoof(x,2))))
        for y in xrange(0,len(numerator_factor_bank)):
            if numerator_factor_bank[y][0] not in possible_primes:
                 if x in possible_factors:
                    possible_factors.remove(x)
    return possible_factors
def reduce_n_times(n,possible_factors, possible_primes, number_of_primes):
    x = 0
    while x < n:
        possible_primes = reduce_possible_primes(possible_factors,number_of_primes)
        possible_factors = reduce_possible_factors(possible_factors,possible_primes)
        x = x+1
    return possible_factors
def create_spoof_matrix(possible_factors,n):
    possible_primes = reduce_possible_primes(possible_factors,n)
    vector = zero_vector(len(possible_primes))
    matrix_list = []
    for x in possible_factors:
        copy_vector = vector[:]
        denominator_list = list(factor(x))
        numerator_list = list(factor(numerator(spoof(x,2))))
        for y in denominator_list:
            copy_vector[possible_primes.index(y[0])] = -2*y[1]
        for y in numerator_list:
            copy_vector[possible_primes.index(y[0])] = y[1]
        matrix_list.append(copy_vector)
    return Matrix(matrix_list)

edit retag flag offensive close merge delete

Comments

Please give us the code, not the picture. Just reedit, copy+paste it into the text, then as long as it is marked press Control+K (or that button with 101 and 010). The following failed for me:

p = MixedIntegerLinearProgram(maximization=True, solver='GLPK' )
x = p.new_variable( integer=True, nonnegative=True )
p.add_constraint( x*spoof_matrix == 0 )

for the obvious reason:

NameError: name 'spoof_matrix' is not defined

dan_fulea ( 2018-03-10 03:52:23 +0100 )edit

I added the code that's running, besides the construction of the matrix which is 56x34, this is literally all that I'm inputting, and I don't know why sage thinks the feasible set is empty since it at least has the trivial solution.

JRHales ( 2018-03-10 04:10:04 +0100 )edit

Could you please provide the code for spoof_matrix ?

tmonteil ( 2018-03-10 12:07:36 +0100 )edit

add a comment

answered 2018-03-12 10:43:50 +0100

dan_fulea
5537 ●5 ●43 ●95

I will try to reproduce the error in a simpler case, with the same error as in the picture.

Let us consider the sample code:

S = matrix( ZZ, 6, 3, [1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 3, 13, 5, 9, 7, 9] )
print S

import traceback

def solveForObjectiveIndex( k ):
    p = MixedIntegerLinearProgram(maximization = True, solver = "GLPK")
    x = p.new_variable(integer = True, nonnegative = True)
    p.add_constraint( x*S == 0)
    p.set_objective( x[k] )
    try:
        return p.solve()
    except Exception:
        traceback.print_exc()
        return

The printed matrix S is:

And now we solve the maximization of $x_k$ problem for possible values of $x$ such that $xS=0$. The kernel of $S$ is:

sage: S.kernel().basis()

[
(3, 13, 0, 5, -1, 0),
(0, 32, 0, 6, -3, 1),
(0, 0, 1, 0, 0, 0)
]

The values of $x$ that can be considered satisfy $x\ge 0$ (on components), if i correctly understand the functionality of the mixed integer linear program framework, and the constraints let us consider only linear combinations with scalars $\ge 0$ because there are unique values $>0$

in column $0$, the entry $3$ at place $(0,0)$,
in column $5$, the entry $1$ at place $(1,5)$,
in column $2$, the entry $1$ at place $(2,2)$).

But the column $4$ tells us that the first two scalars have to be zero. Now we ask for:

sage: solveForObjectiveIndex(0)
0.0
sage: solveForObjectiveIndex(1)
0.0
sage: solveForObjectiveIndex(3)
0.0
sage: solveForObjectiveIndex(4)
0.0
sage: solveForObjectiveIndex(5)
0.0
sage: solveForObjectiveIndex(2)
Traceback (most recent call last):
:::: many lines
MIPSolverException: GLPK: The LP (relaxation) problem has no dual feasible solution

The reason for the maximal value $0$ for indices $k$ among $0,1,3,4,5$ is that the corresponding component is always zero in the set of all $x\ge 0$ with $xS=0$. But for the $2$.nd component we can take arbitrary positive multiples of the basis vector (0, 0, 1, 0, 0, 0), thus the maximal value is $+\infty$. The used algorithm seems to be unprepared for such cases.

Now for the implicitly posted "spoof matrix" with an implicit vector with component $x_{25}=1$ (which does not even appear in the picture) we have the following situation:

spoof_matrix = matrix( ZZ,
[(0, 0, 2, 0, 1, -2, 1, 0, 0, 0, 0, 0, 0, -2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, -2, 0, 0, 0, 0, 0, 1, 0, 0, 0, -2, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-4, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0),
 (0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, -2, 0, 0, 0, 0),
 (0, 0, 1, -2, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 1, 0, 0, 0, 0),
 (1, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 0, -2, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0),
 (-2, -2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 1, 0, 2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 0, 0, -2, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0),
 (1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 2, 0, 0, 0, 1, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-6, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, 2, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 3, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 1, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, -2, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, -2, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, -2, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, -2, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-4, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, -2, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, -4, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 0, 0, 1, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, -2, 0, -2, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 2, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 1, 0, 1, 0, -2, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, -4, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (-4, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 2, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0),
 (-2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0),
 (-4, -2, -2, 0, 2, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 2, 0, 0, -4, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, -2, 0, 1, 0, 1, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (0, 0, 1, 0, 2, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0),
 (1, -2, 2, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
 (1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, -2, 0),
 (0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, -2)] )

(defined over $\mathbb Z$ to have a printable kernel,) and there is exactly one vector in the basis of the kernel with $x_{25}\ne 0$:

sage: S = spoof_matrix
....: for b in S.kernel().basis():
....:     if b[25]:
....:         print b
....: 
(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 28, -32, 4, -30, -11, 16, -16, 0, 30, -22, 0, 0, 0, 6, 0, 0, 0, -24, 3, 16, 0, 0, 0, 0, 0)

The first one is (pythonically) on position $25$. There are some negative values in the list, but one may suppose there is a combination having all components $\ge 0$ and on position $25$ a component $>0$. (This was asked as a bonus.) To obtain a solution, we wil add to the constraints a maximal boundary for $x_{25}$. The code is:

S = spoof_matrix
p = MixedIntegerLinearProgram(maximization = True, solver = "GLPK")
x = p.new_variable(integer = True, nonnegative = True)
p.add_constraint( x*S == 0)
p.add_constraint( x[25] <= 2018 )
p.set_objective( x[25] )
p.solve()
p.get_values(x)

And in the offered long dictionary only the following values of the "corner solution" $x$ are not zero:

 20: 2018.0,
 25: 2018.0,
 29: 4036.0,

(All components are $\ge 0$, the above are mulitples of $1,1,2$. Suitable multiples of this "corner solution" produce bigger and bigger values of the scope function.)

(I hope this answers the question, optimization is not my strength, i feel far away from home.)

edit flag offensive delete link

Comments

Thank you so much! That looks like exactly what I want to do!

JRHales ( 2018-03-14 01:58:53 +0100 )edit

add a comment

Trouble with mixed integer linear programming

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Trouble with mixed integer linear programming edit

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Trouble with mixed integer linear programming