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)
