nikcleju@65: 
nikcleju@65: import numpy
nikcleju@65: import cvxopt
nikcleju@65: import cvxopt.solvers
nikcleju@65: import cvxopt.msk
nikcleju@65: import mosek
nikcleju@65: 
nikcleju@65: def cvxopt_lp(y, A):
nikcleju@65: 
nikcleju@65:     N = A.shape[1]
nikcleju@65:     AA = numpy.hstack((A, -A))
nikcleju@65:     c = numpy.ones((2*N, 1))
nikcleju@65:     
nikcleju@65:     G = numpy.vstack((numpy.zeros((2*N,2*N)), -numpy.eye(2*N)))
nikcleju@65:     h = numpy.zeros((4*N,1))
nikcleju@65: 
nikcleju@65:     # Convert numpy arrays to cvxopt matrices
nikcleju@65:     cvx_c = cvxopt.matrix(c)
nikcleju@65:     cvx_G = cvxopt.matrix(G)
nikcleju@65:     cvx_h = cvxopt.matrix(h)
nikcleju@65:     cvx_AA = cvxopt.matrix(AA)
nikcleju@65:     cvx_y = cvxopt.matrix(y.reshape(y.size,1))
nikcleju@65:     
nikcleju@65:     # Options    
nikcleju@65:     cvxopt.solvers.options['show_progress'] = False
nikcleju@65:     #cvxopt.solvers.options['MOSEK'] = {mosek.iparam.log: 0}
nikcleju@65:     
nikcleju@65:     # Solve
nikcleju@65:     #res = cvxopt.solvers.lp(cvx_c, cvx_G, cvx_h, A=cvx_AA, b=cvx_y, solver='mosek')
nikcleju@65:     res = cvxopt.solvers.lp(cvx_c, cvx_G, cvx_h, A=cvx_AA, b=cvx_y)
nikcleju@65:       
nikcleju@65:     primal = numpy.squeeze(numpy.array(res['x']))
nikcleju@65:     gamma = primal[:N] - primal[N:]
nikcleju@65:     return gamma
nikcleju@65:     
nikcleju@65:     #lb = zeros(2*N,1);
nikcleju@65:     #ub = Inf*ones(2*N,1);
nikcleju@65:     ##[primal,obj,exitflag,output2,dual] = linprog(c,[],[],A,y,lb,ub,[],opt);
nikcleju@65:     ##[primal,~,~,~,~] = linprog(c,[],[],A,aggy,lb,ub);
nikcleju@65:     #[primal,obj,exitflag,output2,dual] = linprog(c,[],[],A,aggy,lb,ub);
nikcleju@65:     #gamma = primal(1:N) - primal((N+1):(2*N));