3
$\begingroup$

I want to optimize the function

$$\min_{X \in \mathbb{S}^{n}_{+}} \mbox{tr} \left( C^T X \right) + \mbox{tr} \left( X^{-1} \right),$$

of which I optimize the equivalent problem

$$\min \mbox{tr}\left(C^T X\right) + \mbox{tr}(Z)$$

$$\text{s.t.} \begin{matrix} \begin{pmatrix} Z & I\\ I & X \\ \end{pmatrix} \end{matrix} \succeq 0 $$

so that it follows the DCP form of cvxpy. Moreover, note that the function has the gradient $C - X^{-2}$. So I define $C$ to be $X_0^{-2}$ for some symmetric positive definite matrix $X_0$. Since the problem is convex, we have that the gradient is $0$ at $X = X_0$. And in this way I can check the performance of the optimization. However, as I put the code in CVXPY, the algorithm suggests that the optimal point is another point. Below is the code and the output.

import cvxpy as cp
import numpy as np
import random
def symm(A):
  #symmetric operator
  return 0.5 *(A + A.T)

np.random.seed(1)
n = 50
CN = 10
data_choice = 'sparse'
solver_choice = 'CG'

D = 1000 * np.diag(np.logspace(-np.log(CN), 0, n))
Q, R = np.linalg.qr(np.random.normal(size = (n,n)))
A = Q @ D @ Q
A = (A.T + A)/2
C= np.linalg.inv(A) @ np.linalg.inv(A)


X = cp.Variable((n, n) , symmetric=True)
Z = cp.Variable((n, n) , symmetric=True)
V = cp.bmat([[Z, np.eye(n)], [np.eye(n), X]])

constraints = [V >> 0]
obj = cp.Minimize(cp.trace(C.T @ X) + cp.trace(Z))
prob = cp.Problem(obj,constraints)

history = prob.solve(verbose = True) 

===============================================================================
                                     CVXPY                                     
                                    v1.1.18                                    
===============================================================================
(CVXPY) Feb 14 12:43:40 AM: Your problem has 5000 variables, 1 constraints, and 0 parameters.
(CVXPY) Feb 14 12:43:40 AM: It is compliant with the following grammars: DCP, DQCP
(CVXPY) Feb 14 12:43:40 AM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Feb 14 12:43:40 AM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
-------------------------------------------------------------------------------
                                  Compilation                                  
-------------------------------------------------------------------------------
(CVXPY) Feb 14 12:43:40 AM: Compiling problem (target solver=SCS).
(CVXPY) Feb 14 12:43:40 AM: Reduction chain: Dcp2Cone -> CvxAttr2Constr -> ConeMatrixStuffing -> SCS
(CVXPY) Feb 14 12:43:40 AM: Applying reduction Dcp2Cone
(CVXPY) Feb 14 12:43:40 AM: Applying reduction CvxAttr2Constr
(CVXPY) Feb 14 12:43:40 AM: Applying reduction ConeMatrixStuffing
(CVXPY) Feb 14 12:43:40 AM: Applying reduction SCS
(CVXPY) Feb 14 12:43:40 AM: Finished problem compilation (took 3.301e-02 seconds).
-------------------------------------------------------------------------------
                                Numerical solver                               
-------------------------------------------------------------------------------
(CVXPY) Feb 14 12:43:40 AM: Invoking solver SCS  to obtain a solution.
------------------------------------------------------------------
           SCS v3.1.0 - Splitting Conic Solver
    (c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 2550, constraints m: 5050
cones:    s: psd vars: 5050, ssize: 1
settings: eps_abs: 1.0e-05, eps_rel: 1.0e-05, eps_infeas: 1.0e-07
      alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
      max_iters: 100000, normalize: 1, warm_start: 0
      acceleration_lookback: 10, acceleration_interval: 10
lin-sys:  sparse-direct
      nnz(A): 2550, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.46e+01  1.00e+00  7.55e+02 -3.51e+02  1.00e-01  3.16e-03 
   250| 1.83e-03  1.19e-03  3.89e-02  3.03e+00  1.00e-01  3.89e-01 
   500| 3.93e-03  1.43e-04  3.54e-03  3.21e+00  4.62e-03  7.70e-01 
   750| 2.25e-03  4.97e-05  2.02e-03  3.00e+00  4.62e-03  1.15e+00 
  1000| 2.02e-03  3.04e-05  1.82e-03  2.99e+00  4.62e-03  1.53e+00 
  1250| 8.64e-04  1.27e-05  7.78e-04  2.90e+00  4.62e-03  1.90e+00 
  1500| 8.21e-04  1.03e-05  7.39e-04  2.90e+00  4.62e-03  2.27e+00 
  1750| 7.89e-04  9.81e-06  7.11e-04  2.90e+00  4.62e-03  2.63e+00 
  2000| 7.63e-04  9.40e-06  6.87e-04  2.90e+00  4.62e-03  3.00e+00 
  2250| 7.42e-04  9.04e-06  6.68e-04  2.89e+00  4.62e-03  3.37e+00 
  2500| 2.22e-04  1.02e-05  2.00e-04  2.86e+00  4.62e-03  3.73e+00 
  2750| 2.44e-04  4.48e-06  2.20e-04  2.86e+00  4.62e-03  4.10e+00 
  3000| 2.34e-04  3.80e-06  2.11e-04  2.86e+00  4.62e-03  4.50e+00 
  3250| 2.30e-04  3.30e-06  2.07e-04  2.86e+00  4.62e-03  4.89e+00 
  3500| 5.20e-03  1.44e-05  1.36e-03  2.84e+00  4.62e-03  5.28e+00 
  3750| 2.24e-04  2.94e-06  2.02e-04  2.86e+00  4.62e-03  5.65e+00 
  4000| 2.21e-04  2.91e-06  1.99e-04  2.86e+00  4.62e-03  6.03e+00 
  4250| 2.19e-04  2.88e-06  1.97e-04  2.86e+00  4.62e-03  6.41e+00 
  4500| 2.16e-04  2.86e-06  1.95e-04  2.86e+00  4.62e-03  6.79e+00 
  4750| 2.14e-04  2.84e-06  1.93e-04  2.86e+00  4.62e-03  7.16e+00 
  5000| 2.12e-04  2.82e-06  1.91e-04  2.86e+00  4.62e-03  7.53e+00 
  5250| 2.10e-04  2.80e-06  1.89e-04  2.86e+00  4.62e-03  7.90e+00 
  5500| 7.40e-05  2.96e-06  6.67e-05  2.85e+00  4.62e-03  8.28e+00 
  5750| 5.86e-05  1.93e-06  5.28e-05  2.85e+00  4.62e-03  8.66e+00 
  6000| 5.93e-05  1.60e-06  5.34e-05  2.85e+00  4.62e-03  9.02e+00 
  6250| 6.10e-05  1.38e-06  5.50e-05  2.85e+00  4.62e-03  9.39e+00 
  6500| 1.43e-04  2.24e-06  7.92e-05  2.85e+00  4.62e-03  9.76e+00 
  6750| 6.32e-05  1.18e-06  5.69e-05  2.85e+00  4.62e-03  1.01e+01 
  7000| 6.37e-05  1.11e-06  5.73e-05  2.85e+00  4.62e-03  1.05e+01 
  7250| 6.38e-05  1.11e-06  5.75e-05  2.85e+00  4.62e-03  1.09e+01 
  7500| 6.38e-05  1.11e-06  5.75e-05  2.85e+00  4.62e-03  1.12e+01 
  7750| 6.37e-05  1.10e-06  5.73e-05  2.85e+00  4.62e-03  1.16e+01 
  8000| 6.34e-05  1.10e-06  5.71e-05  2.85e+00  4.62e-03  1.20e+01 
  8250| 6.31e-05  1.09e-06  5.68e-05  2.85e+00  4.62e-03  1.24e+01 
  8500| 6.27e-05  1.09e-06  5.64e-05  2.85e+00  4.62e-03  1.27e+01 
  8750| 6.22e-05  1.08e-06  5.60e-05  2.85e+00  4.62e-03  1.31e+01 
  9000| 6.18e-05  1.07e-06  5.57e-05  2.85e+00  4.62e-03  1.35e+01 
  9250| 4.37e-03  3.10e-06  7.69e-04  2.84e+00  4.62e-03  1.39e+01 
  9500| 6.09e-05  1.06e-06  5.48e-05  2.85e+00  4.62e-03  1.42e+01 
  9750| 6.04e-05  1.06e-06  5.44e-05  2.85e+00  4.62e-03  1.46e+01 
 10000| 6.00e-05  1.05e-06  5.40e-05  2.85e+00  4.62e-03  1.50e+01 
 10250| 5.95e-05  1.04e-06  5.36e-05  2.85e+00  4.62e-03  1.54e+01 
 10500| 5.91e-05  1.04e-06  5.32e-05  2.85e+00  4.62e-03  1.57e+01 
 10750| 5.87e-05  1.03e-06  5.28e-05  2.85e+00  4.62e-03  1.61e+01 
 11000| 5.82e-05  1.02e-06  5.24e-05  2.85e+00  4.62e-03  1.65e+01 
 11250| 5.78e-05  1.02e-06  5.21e-05  2.85e+00  4.62e-03  1.69e+01 
 11500| 5.74e-05  1.01e-06  5.17e-05  2.85e+00  4.62e-03  1.72e+01 
 11750| 5.70e-05  1.01e-06  5.14e-05  2.85e+00  4.62e-03  1.76e+01 
 12000| 3.48e-03  1.43e-06  6.26e-04  2.84e+00  4.62e-03  1.79e+01 
 12250| 5.63e-05  9.97e-07  5.07e-05  2.85e+00  4.62e-03  1.83e+01 
 12500| 5.60e-05  9.91e-07  5.04e-05  2.85e+00  4.62e-03  1.87e+01 
 12750| 5.56e-05  9.86e-07  5.01e-05  2.85e+00  4.62e-03  1.90e+01 
 13000| 5.53e-05  9.81e-07  4.98e-05  2.85e+00  4.62e-03  1.94e+01 
 13250| 5.50e-05  9.76e-07  4.95e-05  2.85e+00  4.62e-03  1.98e+01 
 13500| 5.47e-05  9.71e-07  4.92e-05  2.85e+00  4.62e-03  2.02e+01 
 13750| 5.44e-05  9.67e-07  4.90e-05  2.85e+00  4.62e-03  2.06e+01 
 13900| 1.01e-05  1.20e-06  9.11e-06  2.84e+00  4.62e-03  2.08e+01 
------------------------------------------------------------------
status:  solved
timings: total: 2.08e+01s = setup: 1.19e-02s + solve: 2.08e+01s
     lin-sys: 9.26e-01s, cones: 1.94e+01s, accel: 1.21e-01s
------------------------------------------------------------------
objective = 2.844050
------------------------------------------------------------------
-------------------------------------------------------------------------------
                                    Summary                                    
-------------------------------------------------------------------------------
(CVXPY) Feb 14 12:44:00 AM: Problem status: optimal
(CVXPY) Feb 14 12:44:00 AM: Optimal value: 2.844e+00
(CVXPY) Feb 14 12:44:00 AM: Compilation took 3.301e-02 seconds
(CVXPY) Feb 14 12:44:00 AM: Solver (including time spent in interface) took 2.081e+01 seconds

I am wondering if anyone can point out if there is any error in the formulation of the problem or the code itself. The optimal value I got from cvxpy is clearly not the optimal value of the function.

np.trace(C.T @ X.value) + np.trace(np.linalg.inv(X.value)) 
2.8435372028588484

np.trace(C.T @ A) + np.trace(np.linalg.inv(A)) 
-0.32275885803468984

So I am wondering if there is any error in my formulation of the problem.

Improve this question
$\endgroup$
4
  • $\begingroup$ @RodrigodeAzevedo It follows from schur complement $\endgroup$
    – The One
    Commented Feb 14, 2022 at 19:00
  • $\begingroup$ From the Schur complement, I get $Z \succeq X^{-1}$, not $Z = X^{-1}$. Could you please give a hint? $\endgroup$
    – Rodrigo de Azevedo
    Commented Feb 14, 2022 at 19:45
  • 1
    $\begingroup$ @Rodrigo de Azevedo The minimization drives that to equality at the optimum. BTW, that is the formulation CVX uses in its function trace_inv, as you can see in lines 21-29 of github.com/cvxr/CVX/blob/master/functions/%40cvx/trace_inv.m . $\endgroup$
    – Mark L. Stone
    Commented Feb 14, 2022 at 20:15
  • $\begingroup$ @MarkL.Stone That makes sense. Thank you. $\endgroup$
    – Rodrigo de Azevedo
    Commented Feb 14, 2022 at 20:33

1 Answer 1

Reset to default
1
$\begingroup$

The "equivalent" problem is equivalent (by Schur Complement of the "equivalent" formulation).

But $A$ is not positive semidefinite. Therefore, $A$ does not satisfy the constraint of the problem. Hence the objective value of $-0.32275885803468984$, obtained from evaluating the objective at $X = A$ does not disprove the optimality of the solution from CVXPY.

Edit: Note, You can't find the optimum by setting the gradient of the objective to zero, as the OP apparently tried, because that neglects the positive semidefiniteness constraint, which contributes to and complicates the optimality condition, which winds up having its own positive semidefiniteness and complementary slackness constraints - hence easier to just solve the original problem numerically.

Improve this answer
$\endgroup$
2
  • $\begingroup$ Perhaps OP meant one of the $Q$s in the definition of $A$ to be $Q^{T}$? Even so it's weird that the entries of $D$ are negative. $\endgroup$
    – Robert Bassett
    Commented Feb 14, 2022 at 19:56
  • $\begingroup$ @Robert Bassett If one of the $Q$ were $Q^T$, and the diagonal of D were nonnegative, then $A$ would be positive semidefinite, and without requiring the symmetrizing operation. But such an $A$ presumably would not have a better objective value than that produced by CVXPY. $\endgroup$
    – Mark L. Stone
    Commented Feb 14, 2022 at 20:19

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.