How to use cvxpy - 10 common examples

lookup = {half_support[i] : i for i in range(len(half_support))}
			constraints = []
			for v,c in coeffs.items():
				if not any(v):
					#constant term gets special treatment
					constraints.append(C[0,0] == coeffs[v] + gamma)
					continue
				#list all (indices of) pairs in half_support, that add up to v
				l = []
				for u in half_support:
					diff = tuple(v[i] - u[i] for i in range(len(v)))
					if diff in half_support:
						l.append((lookup[u],lookup[diff]))
				constraints.append(cvx.Zero(cvx.sum([C[i,j] for i,j in l]) - cvx.expressions.constants.Constant(c)))
			#define the problem
			self.prob_sos_sparse = cvx.Problem(cvx.Minimize(gamma),constraints)
			self.prob_sos = self.prob_sos_sparse
		else:
			n = A.shape[0]
			d = self._degree // 2
			size = binomial(n + 2*d, n)

			#create complete list of monomials, all coefficients initialised with 0
			coeffs = [(list(aux._index_to_vector(i,n,2*d)),0) for i in range(size)]
			#setting the coefficients occurring in A
			for i in range(A.shape[1]):
				index = aux._vector_to_index(A[:,i], 2*d)
				coeffs[index] = (coeffs[index][0],coeffs[index][1] + self.b[i])
			#declare semidefinite matrix C, aim: Z^T * C * Z = p where Z is vector of all monomials
			C = cvx.Variable((binomial(n + d, n), binomial(n + d, n)), PSD = True)

			#construct the constraints

coeffs[tuple(self.A[1:,i])] += self.b[i]
			#create lookup table: vector -> index
			lookup = {half_support[i] : i for i in range(len(half_support))}
			constraints = []
			for v,c in coeffs.items():
				if not any(v):
					#constant term gets special treatment
					constraints.append(C[0,0] == coeffs[v] + gamma)
					continue
				#list all (indices of) pairs in half_support, that add up to v
				l = []
				for u in half_support:
					diff = tuple(v[i] - u[i] for i in range(len(v)))
					if diff in half_support:
						l.append((lookup[u],lookup[diff]))
				constraints.append(cvx.Zero(cvx.sum([C[i,j] for i,j in l]) - cvx.expressions.constants.Constant(c)))
			#define the problem
			self.prob_sos_sparse = cvx.Problem(cvx.Minimize(gamma),constraints)
			self.prob_sos = self.prob_sos_sparse
		else:
			n = A.shape[0]
			d = self._degree // 2
			size = binomial(n + 2*d, n)

			#create complete list of monomials, all coefficients initialised with 0
			coeffs = [(list(aux._index_to_vector(i,n,2*d)),0) for i in range(size)]
			#setting the coefficients occurring in A
			for i in range(A.shape[1]):
				index = aux._vector_to_index(A[:,i], 2*d)
				coeffs[index] = (coeffs[index][0],coeffs[index][1] + self.b[i])
			#declare semidefinite matrix C, aim: Z^T * C * Z = p where Z is vector of all monomials
			C = cvx.Variable((binomial(n + d, n), binomial(n + d, n)), PSD = True)

def test_basic_gp(self):
        tf.random.set_seed(243)

        x = cp.Variable(pos=True)
        y = cp.Variable(pos=True)
        z = cp.Variable(pos=True)

        a = cp.Parameter(pos=True, value=2.0)
        b = cp.Parameter(pos=True, value=1.0)
        c = cp.Parameter(value=0.5)

        objective_fn = 1/(x*y*z)
        constraints = [a*(x*y + x*z + y*z) <= b, x >= y**c]
        problem = cp.Problem(cp.Minimize(objective_fn), constraints)
        problem.solve(cp.SCS, gp=True, eps=1e-12)

        layer = CvxpyLayer(
            problem, parameters=[a, b, c], variables=[x, y, z], gp=True)
        a_tf = tf.Variable(2.0, dtype=tf.float64)
        b_tf = tf.Variable(1.0, dtype=tf.float64)

"""Solve a group elastic net problem with cvx.

    Arguments:
        As: list of tensors.
        y: tensor.
        l_1, l_2: floats.
        ns: iterable.
    """
    # Convert everything to numpy
    dtype = As[0].dtype
    As = [A_j.cpu().numpy() for A_j in As]
    y = y.cpu().numpy()
    ns = np.array([int(n) for n in ns])

    # Create the b variables.
    b_0 = cvx.Variable()
    bs = []
    for _, n_j in zip(As, ns):
        bs.append(cvx.Variable(n_j))

    # Form g(b).
    Ab = sum(A_j * b_j for A_j, b_j in zip(As, bs))
    m = As[0].shape[0]
    g_b = cvx.square(cvx.norm(y - b_0 - Ab)) / (2 * m)

    # Form h(b).
    h_b = sum(
        np.sqrt(n_j) * (l_1 * cvx.norm(b_j) + l_2 * cvx.square(cvx.norm(b_j)))
        for n_j, b_j in zip(ns, bs)
    )

    # Build the optimization problem.

#    np.log(trAtA)
                       ])

    # Add row to matrix A
    A = spa.vstack((A, spa.csc_matrix(A_temp)), 'csc')

    # Add bounds on v
    l = np.log(problem['rho_min'].iloc[0])
    u = np.log(problem['rho_max'].iloc[0])
    v_l = np.append(v_l, l)
    v_u = np.append(v_u, u)

#
# Define CVXPY problem
alpha = cvxpy.Variable(n_params)
v = cvxpy.Variable(n_problems)

constraints = [v_l <= v, v <= v_u]
cost = cvxpy.norm(A * alpha - v)

objective = cvxpy.Minimize(cost)

problem = cvxpy.Problem(objective, constraints)

# Solve problem
print("Solving problem with CVXPY and GUROBI")
problem.solve(solver=cvxpy.MOSEK, verbose=True)
print("Solution status: %s" % problem.status)

# Get learned alpha
alpha_fit = np.asarray(alpha.value).flatten()

device = A_j.device
    dtype = A_j.dtype
    r_j = r_j.cpu().numpy()
    A_j = A_j.cpu().numpy()

    # Create the b_j variable.
    b_j = cvx.Variable(A_j.shape[1])

    # Form the objective.
    q_j = r_j - A_j * b_j
    obj_fun = cvx.square(cvx.norm(q_j)) / (2.0 * m)
    obj_fun += a_1_j * cvx.norm(b_j) + (a_2_j / 2.0) * cvx.square(cvx.norm(b_j))

    # Build the optimization problem.
    obj = cvx.Minimize(obj_fun)
    problem = cvx.Problem(obj, constraints=None)

    problem.solve(solver="CVXOPT", verbose=False)
    b_j = np.asarray(b_j.value)
    return torch.from_numpy(b_j).to(device, dtype)

"""
        import numpy

        # Problem data.
        n = 15
        m = 10
        numpy.random.seed(1)
        A = numpy.random.randn(n, m)
        b = numpy.random.randn(n)
        # gamma must be positive due to DCP rules.
        gamma = cp.Parameter(nonneg=True)

        # Construct the problem.
        x = cp.Variable(m)
        error = cp.sum_squares(A*x - b)
        obj = cp.Minimize(error + gamma*cp.norm(x, 1))
        prob = cp.Problem(obj)

        # Construct a trade-off curve of ||Ax-b||^2 vs. ||x||_1
        sq_penalty = []
        l1_penalty = []
        x_values = []
        gamma_vals = numpy.logspace(-4, 6, 10)

        start = time.time()
        for val in gamma_vals:
            gamma.value = val
            prob.solve(solver=cp.SCS, warm_start=True, use_indirect=True)
            # Use expr.value to get the numerical value of
            # an expression in the problem.
            sq_penalty.append(error.value)
            l1_penalty.append(cp.norm(x, 1).value)

def test_too_many_variables(self):
        x = cp.Variable(1)
        y = cp.Variable(1)
        lam = cp.Parameter(1, nonneg=True)
        objective = lam * cp.norm(x, 1)
        prob = cp.Problem(cp.Minimize(objective))
        with self.assertRaisesRegex(ValueError, 'Argument `variables`.*'):
            CvxpyLayer(prob, [lam], [x, y])  # noqa: F841

def test_warmstart(self):
        """Testing warmstart LASSO with SCS.
        """
        import numpy

        # Problem data.
        n = 15
        m = 10
        numpy.random.seed(1)
        A = numpy.random.randn(n, m)
        b = numpy.random.randn(n)
        # gamma must be positive due to DCP rules.
        gamma = cp.Parameter(nonneg=True)

        # Construct the problem.
        x = cp.Variable(m)
        error = cp.sum_squares(A*x - b)
        obj = cp.Minimize(error + gamma*cp.norm(x, 1))
        prob = cp.Problem(obj)

        # Construct a trade-off curve of ||Ax-b||^2 vs. ||x||_1
        sq_penalty = []
        l1_penalty = []
        x_values = []
        gamma_vals = numpy.logspace(-4, 6, 10)

        start = time.time()
        for val in gamma_vals:
            gamma.value = val

np.random.seed(243)
        N, n = 10, 2

        def sigmoid(z):
            return 1 / (1 + np.exp(-z))

        X_np = np.random.randn(N, n)
        a_true = np.random.randn(n, 1)
        y_np = np.round(sigmoid(X_np @ a_true + np.random.randn(N, 1) * 0.5))

        X_tf = tf.Variable(X_np)
        lam_tf = tf.Variable(1.0 * tf.ones(1))

        a = cp.Variable((n, 1))
        X = cp.Parameter((N, n))
        lam = cp.Parameter(1, nonneg=True)
        y = y_np

        log_likelihood = cp.sum(
            cp.multiply(y, X @ a) -
            cp.log_sum_exp(cp.hstack([np.zeros((N, 1)), X @ a]).T, axis=0,
                           keepdims=True).T
        )
        prob = cp.Problem(
            cp.Minimize(-log_likelihood + lam * cp.sum_squares(a)))
        fit_logreg = CvxpyLayer(prob, [X, lam], [a])

        with tf.GradientTape(persistent=True) as tape:
            weights = fit_logreg(X_tf, lam_tf, solver_args={'eps': 1e-8})[0]
            summed = tf.math.reduce_sum(weights)
        grad_X_tf, grad_lam_tf = tape.gradient(summed, [X_tf, lam_tf])