Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def compare_opt(self, a_ub, b_ub, c):
'compare cvx opt versus our glpk interface'
# make sure we're using floats not ints
a_ub = [[float(x) for x in row] for row in a_ub]
b_ub = [float(x) for x in b_ub]
c = [float(x) for x in c]
num_vars = len(a_ub[0])
# solve it with cvxopt
options = {'show_progress': False}
sol = cvxopt.solvers.lp(cvxopt.matrix(c), cvxopt.matrix(a_ub).T, cvxopt.matrix(b_ub), options=options)
if sol['status'] != 'optimal':
raise RuntimeError("cvxopt LP failed: {}".format(sol['status']))
res_cvxopt = [float(n) for n in sol['x']]
# solve it with the glpk <-> hylaa interface
lp = LpInstance(num_vars, num_vars)
lp.set_init_constraints(csr_matrix(np.array(a_ub, dtype=float)), np.array(b_ub, dtype=float))
lp.set_no_output_constraints()
lp.update_basis_matrix(np.identity(num_vars))
res_glpk = np.zeros(num_vars)
lp.minimize(np.array(c, dtype=float), res_glpk)
h.append(p_interval.high)
h = matrix(h)
#print h
# A matrix, to be used as A:
# form: same as c
A = matrix(1.0, (1,n))
#print A
# b matrix
# form: 1x1 with entry 1
b = matrix(1.0, (1,1))
#print b
solvers.options["show_progress"] = False
P = solvers.lp(V, G, h, A, b)['x']
#for p in P:
# print p
#print
return P
c = cvxopt.matrix(coef, tc='d')
# for positivity inequalities
G = cvxopt.spdiag(cvxopt.matrix(-np.ones(n_variables)))
#G = cvxopt.matrix(-np.eye(n_variables))
h = cvxopt.matrix(np.zeros(n_variables)) # for positivity inequalities
# unary and pairwise summation constratints
A = cvxopt.spmatrix(data, I, J)
assert(n_constraints == A.size[0])
b_ = np.zeros(A.size[0]) # zeros for pairwise summation constraints
b_[:n_nodes] = 1 # ones for unary summation constraints
b = cvxopt.matrix(b_)
# don't be verbose.
show_progress_backup = cvxopt.solvers.options.get('show_progress', False)
cvxopt.solvers.options['show_progress'] = False
result = cvxopt.solvers.lp(c, G, h, A, b)
cvxopt.solvers.options['show_progress'] = show_progress_backup
x = np.array(result['x'])
unary_variables = x[:n_nodes * n_states].reshape(n_nodes, n_states)
pairwise_variables = x[n_nodes * n_states:].reshape(n_edges, n_states ** 2)
assert((np.abs(unary_variables.sum(axis=1) - 1) < 1e-4).all())
assert((np.abs(pairwise_variables.sum(axis=1) - 1) < 1e-4).all())
return unary_variables, pairwise_variables, result['primal objective']
# try:
# import mosek
# cvxopt.solvers.options['MOSEK'] = {mosek.iparam.log: 0}
# except:
# warnings.warn('No Mosek import possible!')
# Do the optimization
# (-c_i x_i,c_i x_i, -c_i,-e_i) * (w+,w-,b,t) \leq -1
# (c_i x_i,-c_i x_i, c_i,-e_i) * (w+,w-,b,t) \leq R
# (-w+,-w-,-t) \leq 0
self._log("Classifier under construction")
if not self.max_time is None:
cvxopt.solvers.options['LPX_K_TMLIM'] = int(self.max_time)
self.sol = cvxopt.solvers.lp(
c,
cvxopt.matrix(self.create_problem_matrix(n, e)), h, solver='glpk')
self._log("Construction complete")
model = []
self.calculate_classification_vector(model)
if self.debug:
self.calculate_slack_variables(model)
inner_dual_sol = self.sol['z'][:n]
dual_inner_weights = numpy.multiply(numpy.array(inner_dual_sol).T,
numpy.array(self.bi))[0]
self.max_inner_weight = []
self.max_inner_weight.append(-dual_inner_weights.min())
self.max_inner_weight.append(dual_inner_weights.max())
`C, D, b`: Polytope parameters
Output:
`x_0, y_0`: The chebyshev centra
`boolean`: True if a point could be found, False otherwise'''
d = C.shape[1]
k = D.shape[1]
A = np.hstack([C,D])
dim = np.shape(A)[1]
c = -np.r_[np.zeros(dim),1]
norm2 = np.sqrt(np.sum(A*A, axis=1))
G = np.c_[A, norm2]
solvers.options['show_progress']=False
solvers.options['LPX_K_MSGLEV'] = 0
sol = solvers.lp(matrix(c), matrix(G), matrix(b), None, None, lp_solver)
if sol['status'] == "optimal":
opt = np.array(sol['x'][0:-1]).flatten()
return opt[range(0,d)], opt[range(d,d+k)], True
else:
return np.zeros(d), np.zeros(k), False
solver : string, optional
Backend LP solver to call.
Returns
-------
succ : bool
Success boolean.
z : (3,) array, or 0
Maximum vertex of the polygon in the direction `vdir`, or 0 in case of
solver failure.
"""
lp_q, lp_Gextended, lp_hextended, lp_A, lp_b = lp
lp_q[-2] = -vdir[0]
lp_q[-1] = -vdir[1]
try:
sol = cvxopt.solvers.lp(
lp_q, lp_Gextended, lp_hextended, lp_A, lp_b, solver=solver)
if sol['status'] == 'optimal':
z = sol['x']
z = array(z).reshape((lp_q.size[0], ))
return True, z[-2:]
else:
warn("Failed with status %s\n" % sol['status'])
return False, 0
except Exception as inst:
warn("Catched exception: {}".format(inst))
return False, 0
sage: sol['x']
(45.0..., 6.25, 1.0...)
"""
from cvxopt.base import matrix as m
from cvxopt import solvers
solvers.options['show_progress']=False
if solver=='glpk':
from cvxopt import glpk
glpk.options['LPX_K_MSGLEV'] = 0
c_=m(c.base_extend(RDF).numpy())
G_=m(G.base_extend(RDF).numpy())
h_=m(h.base_extend(RDF).numpy())
if A!=None and b!=None:
A_=m(A.base_extend(RDF).numpy())
b_=m(b.base_extend(RDF).numpy())
sol=solvers.lp(c_,G_,h_,A_,b_,solver=solver)
else:
sol=solvers.lp(c_,G_,h_,solver=solver)
status=sol['status']
# CHANGED:
#if status != 'optimal':
if status not in ('optimal', 'dual infeasible'):
return {'x':None,'s':None,
#'y':None, 'z':None,
'status':status}
x=vector(RDF,list(sol['x']))
s=vector(RDF,list(sol['s']))
#y=vector(RDF,list(sol['y']))
#z=vector(RDF,list(sol['z']))
#return {'primal objective':sol['primal objective'],'x':x,'s':s,'y':y,
# 'z':z,'status':status}
return {'x': x, 's': s, 'status': status}