How to use the cvxpy.square function in cvxpy

np.random.seed(1)
HR = np.random.randn(m, n)
HI = np.random.randn(m, n)
A = np.hstack((HR, HI))
B = np.hstack((-HI, HR))

GR = np.random.randn(l, n)
GI = np.random.randn(l, n)
C = np.hstack((GR, GI))
D = np.hstack((-GI, GR))

x = cvx.Variable(2*n)
obj = cvx.Minimize(cvx.sum_squares(x))
cons = [
    cvx.square(A*x) + cvx.square(B*x) >= tau,
    cvx.square(C*x) + cvx.square(D*x) <= eta
]
prob = cvx.Problem(obj, cons)
qcqp = QCQP(prob)

# sample from the semidefinite relaxation
qcqp.suggest(SDR)
print("SDR-based lower bound: %.3f" % qcqp.sdr_bound)

# f_cd, v_cd = qcqp.improve(COORD_DESCENT)
# print("Coordinate descent: objective %.3f, violation %.3f" % (f_cd, v_cd))
# f_cd, v_cd = qcqp.improve(ADMM, rho=np.sqrt(m+l), phase1=False)
# print("Coordinate descent: objective %.3f, violation %.3f" % (f_cd, v_cd))

# # SDR solution is cached and not solved again
# qcqp.suggest(SDR)

cvxpy.reshape(self.par['A_bar'][:, k], (m.n_x, m.n_x)) *
            self.var['X'][:, k] +
            cvxpy.reshape(self.par['B_bar'][:, k], (m.n_x, m.n_u)) *
            self.var['U'][:, k] +
            cvxpy.reshape(self.par['C_bar'][:, k], (m.n_x, m.n_u)) *
            self.var['U'][:, k + 1] +
            self.par['S_bar'][:, k] * self.var['sigma'] +
            self.par['z_bar'][:, k] +
            self.var['nu'][:, k]
            for k in range(K - 1)
        ]

        # Trust regions:
        dx = cvxpy.sum(cvxpy.square(
            self.var['X'] - self.par['X_last']), axis=0)
        du = cvxpy.sum(cvxpy.square(
            self.var['U'] - self.par['U_last']), axis=0)
        ds = self.var['sigma'] - self.par['sigma_last']
        constraints += [cvxpy.norm(dx + du, 1) <= self.var['delta_norm']]
        constraints += [cvxpy.norm(ds, 'inf') <= self.var['sigma_norm']]

        # Flight time positive:
        constraints += [self.var['sigma'] >= 0.1]

        # Objective:
        sc_objective = cvxpy.Minimize(
            self.par['weight_sigma'] * self.var['sigma'] +
            self.par['weight_nu'] * cvxpy.norm(self.var['nu'], 'inf') +
            self.par['weight_delta'] * self.var['delta_norm'] +
            self.par['weight_delta_sigma'] * self.var['sigma_norm']
        )

.. math::

            minimize(norm((signal)- (phi*fe)))
        """

        # Create four scalar optimization variables.
        fe = cvx.Variable(2)
        # Create four constraints.
        constraints = [cvx.sum(fe) == 1,
                       fe[0] >= 0.011,
                       fe[1] >= 0.011,
                       fe[0] <= 0.29,
                       fe[1] <= 0.89]

        # Form objective.
        obj = cvx.Minimize(cvx.sum(cvx.square(phi * fe - signal)))

        # Form and solve problem.
        prob = cvx.Problem(obj, constraints)
        prob.solve()  # Returns the optimal value.
        return np.array(fe.value)

"""

    warnings.warn('This objective is will try to maximize a convex function.')

    obj = 0

    c_vars = model_vars['commenters']
    v_vars = model_vars['voters']
    for topic_id, p_id, c_id, v_id, p_v, c_v in \
            df[['ArticleTopic', 'ParentCommenterId', 'ChildCommenterId',
                'VoterId', 'ParentVote', 'ChildVote']].values:
        if p_v < 0:
            obj += CVX.square(c_vars[commenter_id_to_idx[p_id]][topic_id_to_idx[topic_id]] -
                              v_vars[voter_id_to_idx[v_id]][topic_id_to_idx[topic_id]])
        if c_v < 0:
            obj += CVX.square(c_vars[commenter_id_to_idx[c_id]][topic_id_to_idx[topic_id]] -
                              v_vars[voter_id_to_idx[v_id]][topic_id_to_idx[topic_id]])

    return CVX.Maximize(obj)

import cvxpy as cvx
import dmcp

n = 10
m = 15
I = np.eye(n)
C = np.random.randn(m,n) # channel matrix
U, S, V = np.linalg.svd(C)
U = U[:,0:n]
sigma_e = cvx.Parameter(nonneg=True) # additive Gaussian noise

A = cvx.Variable((n,n)) # pre-coder
B = cvx.Variable((n,m)) # equalizer

#cost = square(norm(B*C*A-I,'fro'))/2+square(sigma_e)*square(norm(B,'fro')) # P(xi=0) = P(xi=1) = 0.5
cost = cvx.sum(cvx.square(B*C*A-I))/2 + cvx.square(sigma_e)*cvx.sum(cvx.square(B))
obj = cvx.Minimize(cost/(m*n))
prob = cvx.Problem(obj, [cvx.norm(A,'fro')<=10])

SNR = np.power(10,np.linspace(-2,1,20))
sigma_e_value = np.sqrt(float(n)/2/m/SNR)

MSE = []
MSE_lin = []
MSE_noprox = []
BER = []
BER_lin = []
BER_noprox = []
for value in sigma_e_value:
    sigma_e.value = value
    B.value = np.dot(np.dot(np.transpose(V), np.diag(float(1)/np.sqrt(S))), np.transpose(U))
    A.value = np.dot(np.dot(np.transpose(V), np.diag(float(1)/np.sqrt(S))), V)

:param var_slack: list of slack variables
    :param lambd: proximal operator parameter
    :return: a problem with proximal operator
    """
    new_cost = prob.objective.expr
    new_constr = prob.constraints
    slack_id = [var.id for var in var_slack]

    #Add proximal variables
    prob_variables = prob.variables()
    prob_variables.sort(key = lambda x:x.id)
    for var in prob_variables:
        # add quadratic terms for all variables that are not slacks
        if not var.id in slack_id:
            if prob.objective.NAME == 'minimize':
                new_cost = new_cost + cvx.square(cvx.norm(var - var.value,'fro'))/2/lambd
            else:
                new_cost = new_cost - cvx.square(cvx.norm(var - var.value,'fro'))/2/lambd

    # Define proximal problem
    if prob.objective.NAME == 'minimize':
        new_prob = cvx.Problem(cvx.Minimize(new_cost), new_constr)
    else: # maximize
        new_prob = cvx.Problem(cvx.Maximize(new_cost), new_constr)
    return new_prob

n = 25
np.random.seed(1)

# Make adjacency matrix.
p = 0.2
W = np.asmatrix(np.random.uniform(low=0.0, high=1.0, size=(n, n)))
for i in range(n):
    W[i, i] = 1
    for j in range(i+1, n):
        W[j, i] = W[i, j]
W = (W < p).astype(float)

x = cvx.Variable(n)
obj = 0.25*(cvx.sum_entries(W) - cvx.quad_form(x, W))
cons = [cvx.square(x) == 1]
prob = cvx.Problem(cvx.Maximize(obj), cons)

# Objective function
f = lambda x: 0.25*(np.sum(W) - (x_round.T*W*x_round)[0, 0])
def violation(x):
    return np.max(np.abs(np.square(x) - 1))

# SDP-based upper bound
ub = prob.solve(method='sdp-relax', solver=cvx.MOSEK)
print ('Upper bound: %.3f' % ub)

# Lower bounds
print ('(objective, maximum violation):')
f_dccp = prob.solve(method='qcqp-dccp', use_sdp=False, solver=cvx.MOSEK, num_samples=10, tau=1)
print ('  Convex-concave programming: (%.3f, %.3f)' % (f_dccp, violation(x.value)))
x_round = np.sign(np.random.randn(n, 1))