It seems that square_pos() is wrong

Here is my code.

cvx_begin
cvx_solver mosek_2
%优化变量, x1,y1为新的轨迹
variables t x1(101) y1(101)
%设置变量关系
maximize(t)
subject to
e1 = Puh0./a1;
e2 = Pu
h0./a2;
e3 = Pu*h0./a3;
c1 = P1.*h0./b1;
c2 = P2.*h0./b2;
c3 = P3.h0./b3;
%各用户可达速率的定义
R1ulb = a1.log(1+e1./(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2))./log(2)-a1./log(2).e1./((H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2).(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2+e1)).(square_pos(norms([x1 y1]-repmat([s1x s1y],101,1),[],2))-norms([x y]-repmat([s1x s1y],101,1),[],2).^2);
R2ulb = a2.log(1+e2./(H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2))./log(2)-a2./log(2).e2./((H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2).(H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2+e2)).(square_pos(norms([x1 y1]-repmat([s2x s2y],101,1),[],2))-norms([x y]-repmat([s2x s2y],101,1),[],2).^2);
R3ulb = a3.log(1+e3./(H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2))./log(2)-a3./log(2).e3./((H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2).(H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2+e3)).(square_pos(norms([x1 y1]-repmat([s3x s3y],101,1),[],2))-norms([x y]-repmat([s3x s3y],101,1),[],2).^2);
R1vlb = b1.log(1+c1./(H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2))./log(2)-b1./log(2).c1./((H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2).(H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2+c1)).(square_pos(norms([x1 y1]-repmat([d1x d1y],101,1),[],2))-norms([x y]-repmat([d1x d1y],101,1),[],2).^2);
R2vlb = b2.log(1+c2./(H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2))./log(2)-b2./log(2).c2./((H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2).(H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2+c2)).(square_pos(norms([x1 y1]-repmat([d2x d2y],101,1),[],2))-norms([x y]-repmat([d2x d2y],101,1),[],2).^2);
R3vlb = b3.log(1+c3./(H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2))./log(2)-b3./log(2).c3./((H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2).(H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2+c3)).(square_pos(norms([x1 y1]-repmat([d3x d3y],101,1),[],2))-norms([x y]-repmat([d3x d3y],101,1),[],2).^2);
%约束条件(22a)和(22b),因为有六个用户,所以分为六个式子0.
B0/(101
Rt)sum(R1ulb) >= t;
B0/(101
Rt)sum(R2ulb) >= t;
B0/(101
Rt)sum(R3ulb) >= t;
B0/(101
Rt)sum(R1vlb) >= t;
B0/(101
Rt)sum(R2vlb) >= t;
B0/(101
Rt)*sum(R3vlb) >= t;
%约束条件(18h)和(18i)
for i=1:100
norm([x1(i+1)-x1(i),y1(i+1)-y1(i)]) <= Dmax
end
x1(101) == x1(1);
y1(101) == y1(1);

cvx_end

After I got the result, x1, y1 is solved. Then I replace square_pos() in R1ulb with .^2. The result should be the same, but it wasn’t. I am very confused. The norm is always positive so square_pos() should equal to .^2. But it wasn’t. Can someone help me?
%E6%89%B9%E6%B3%A8%202020-03-06%20134202
The orange line is R1ulb with square_pos() in cvx. The blue line is R1ulb after i replace square_pos() with .^2.

Can you reduce this to a simpler problem which still exhibits the phenomenon, and show the exact code used for both versions/

you can focus on R1ulb. In my problem, the variables in cvx are x1 and y1, both of them are 101 demention vectors. To make it easier to understand, R1ulb can be simplyfied as R1ulb=M-N.(square_pos(norms([x1 y1]-repmat([s1x s1y],101,1),[],2))-K), where M,N,K(M,N,K are also 101 demention vectors) are given. So I run the code, I got the result t,x1,y1. But the result doesn’t seem right. Then I write Riulb_new=M-N.(norms([x1 y1]-repmat([s1x s1y],101,1),[],2).^2-K), which simply replace square_pos() with .^2. I think R1ulb and R1ulb_new should be the same, because square_pos() is just cvx version of .^2. But R1ulb and R1ulb_new are so different, just like the graph shows. I wonder where went wrong, if this is a bug in cvx?

you can focus on R1ulb. In my problem, the variables in cvx are x1 and y1, both of them are 101 demention vectors. To make it easier to understand, R1ulb can be simplyfied as R1ulb=M-N. *(square_pos(norms([x1 y1]-repmat([s1x s1y],101,1),[],2))-K), where M,N,K(M,N,K are also 101 demention vectors) are given. So I run the code, I got the result t,x1,y1. But the result doesn’t seem right. Then I write Riulb_new=M-N. *(norms([x1 y1]-repmat([s1x s1y],101,1),[],2).^2-K), which simply replace square_pos() with .^2. I think R1ulb and R1ulb_new should be the same, because square_pos() is just cvx version of .^2. But R1ulb and R1ulb_new are so different, just like the graph shows. I wonder where went wrong, if this is a bug in cvx?

What is the output of cvx_version ?

It seems that you must be usiing CVX 3. beta. otherwise, norm(...)^2 would not be accepted by CVX. Unfortunately, CVX 3.0beta is riddled with bugs, so I recommend you install CVX 2.2, in which case you will have to use square_pos, because norm(...)^2 is not allowed in CVX 2.2…

My cvx version is 2.1. Riulb_new=M-N.*(norms([x1 y1]-repmat([s1x s1y],101,1),[],2).^2-K) is out of the cvx area just to check square_pos() in R1ulb is right or not, because after I got the result , x1 and y1 can be assumed as given.

You are making it difficult for me to see clearly what you have done. That is why I asked for complete code for both versions.

In any event, if you are using CVX 2.1, the norm(...)^2 must have been applied after CVX ended. In that case, presuming CVX Status is Solved, you can only be guaranteed that CVX variables are at their optimal values. Expressions computed in terms of CVX variables may not be at their optimal values after CVX ends, even though the optimization problem was solved correctly. So to compute expressions after CVX ends, you must doing ti starting only from CVX variables, not expressions, unless they are in turn computed after CVX ends using only CVX variables.

But I cannot answer your question, because you have not shown clearly what you have done.

And please use the Preformatted text icon (6th from the left) so that your code displays properly.

cvx_begin
cvx_solver mosek_2
variables t x1(101) y1(101)
maximize(t)
subject to
    e1 = Pu*h0./a1;
    e2 = Pu*h0./a2;
    e3 = Pu*h0./a3;
    c1 = P1.*h0./b1;
    c2 = P2.*h0./b2;
    c3 = P3.*h0./b3;
    R1ulb = a1.*log(1+e1./(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2))./log(2)-a1./log(2).*e1./((H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2+e1)).*(square_pos(norms([x1 y1]-repmat([s1x s1y],101,1),[],2))-norms([x y]-repmat([s1x s1y],101,1),[],2).^2);
    R2ulb = a2.*log(1+e2./(H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2))./log(2)-a2./log(2).*e2./((H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([s2x s2y],101,1),[],2).^2+e2)).*(square_pos(norms([x1 y1]-repmat([s2x s2y],101,1),[],2))-norms([x y]-repmat([s2x s2y],101,1),[],2).^2);
    R3ulb = a3.*log(1+e3./(H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2))./log(2)-a3./log(2).*e3./((H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([s3x s3y],101,1),[],2).^2+e3)).*(square_pos(norms([x1 y1]-repmat([s3x s3y],101,1),[],2))-norms([x y]-repmat([s3x s3y],101,1),[],2).^2);
    R1vlb = b1.*log(1+c1./(H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2))./log(2)-b1./log(2).*c1./((H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([d1x d1y],101,1),[],2).^2+c1)).*(square_pos(norms([x1 y1]-repmat([d1x d1y],101,1),[],2))-norms([x y]-repmat([d1x d1y],101,1),[],2).^2);
    R2vlb = b2.*log(1+c2./(H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2))./log(2)-b2./log(2).*c2./((H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([d2x d2y],101,1),[],2).^2+c2)).*(square_pos(norms([x1 y1]-repmat([d2x d2y],101,1),[],2))-norms([x y]-repmat([d2x d2y],101,1),[],2).^2);
    R3vlb = b3.*log(1+c3./(H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2))./log(2)-b3./log(2).*c3./((H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([d3x d3y],101,1),[],2).^2+c3)).*(square_pos(norms([x1 y1]-repmat([d3x d3y],101,1),[],2))-norms([x y]-repmat([d3x d3y],101,1),[],2).^2);
    B0/(101*Rt)*sum(R1ulb) >= t;
    B0/(101*Rt)*sum(R2ulb) >= t;
    B0/(101*Rt)*sum(R3ulb) >= t;
    B0/(101*Rt)*sum(R1vlb) >= t;
    B0/(101*Rt)*sum(R2vlb) >= t;
    B0/(101*Rt)*sum(R3vlb) >= t;
    for i=1:100
        norm([x1(i+1)-x1(i),y1(i+1)-y1(i)]) <= Dmax
    end
    x1(101) == x1(1);
    y1(101) == y1(1);
cvx_end

R1ulb_2 = a1.*log(1+e1./(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2))./log(2)-a1./log(2).*e1./((H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2).*(H^2+norms([x y]-repmat([s1x s1y],101,1),[],2).^2+e1)).*(norms([x1 y1]-repmat([s1x s1y],101,1),[],2).^2-norms([x y]-repmat([s1x s1y],101,1),[],2).^2);

My complete code is here. What confused me is that R1ulb_2 should be close to R1ulb after CVX process. But they are so different as the previous graph shows. The status is Inaccurate/Solved.

Please show full log output from Mosek.

R1ulb is a CVX expression involving the variables x and y1. Therefore, its value after CVX completes may not be correct (even though the optimization problem is solved correctly (leaving aside “inaccurate” status)). R1ulb_2 should be the correct value (at optimum) for R1ulb, because expressions must be recomputed, starting from CVX variables, after CVX completion, ir order to be guaranteed to have the correct value.

That resolves your original concern. But please go ahead and provided the Mosek log output so that @Michal_Adamaszek can assess the “inaccurate” status of the solution.

Calling Mosek_2 9.1.10: 5255 variables, 2530 equality constraints
------------------------------------------------------------

MOSEK Version 9.1.10 (Build date: 2020-1-1 23:02:09)
Copyright (c) MOSEK ApS, Denmark. WWW: mosek.com
Platform: Windows/64-X86

MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (304) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (307) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (310) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (313) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (316) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (319) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (322) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (325) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (328) of matrix 'A'.
MOSEK warning 710: #1 (nearly) zero elements are specified in sparse col '' (331) of matrix 'A'.
Warning number 710 is disabled.
Problem
  Name                   :                 
  Objective sense        : min             
  Type                   : CONIC (conic optimization problem)
  Constraints            : 2530            
  Cones                  : 1312            
  Scalar variables       : 5255            
  Matrix variables       : 0               
  Integer variables      : 0               

Optimizer started.
Presolve started.
Linear dependency checker started.
Linear dependency checker terminated.
Eliminator started.
Freed constraints in eliminator : 100
Eliminator terminated.
Eliminator started.
Freed constraints in eliminator : 118
Eliminator terminated.
Eliminator - tries                  : 2                 time                   : 0.00            
Lin. dep.  - tries                  : 1                 time                   : 0.00            
Lin. dep.  - number                 : 0               
Presolve terminated. Time: 0.03    
Problem
  Name                   :                 
  Objective sense        : min             
  Type                   : CONIC (conic optimization problem)
  Constraints            : 2530            
  Cones                  : 1312            
  Scalar variables       : 5255            
  Matrix variables       : 0               
  Integer variables      : 0               

Optimizer  - threads                : 4               
Optimizer  - solved problem         : the primal      
Optimizer  - Constraints            : 969
Optimizer  - Cones                  : 743
Optimizer  - Scalar variables       : 2555              conic                  : 2228            
Optimizer  - Semi-definite variables: 0                 scalarized             : 0               
Factor     - setup time             : 0.00              dense det. time        : 0.00            
Factor     - ML order time          : 0.00              GP order time          : 0.00            
Factor     - nonzeros before factor : 5690              after factor           : 8709            
Factor     - dense dim.             : 0                 flops                  : 1.13e+05        
ITE PFEAS    DFEAS    GFEAS    PRSTATUS   POBJ              DOBJ              MU       TIME  
0   4.1e+03  1.0e+00  1.0e+04  0.00e+00   0.000000000e+00   -1.016050000e+04  1.0e+00  0.05  
1   6.9e+02  1.7e-01  4.2e+03  -1.00e+00  4.344302542e+00   -1.015079909e+04  1.7e-01  0.13  
2   1.6e+02  4.0e-02  2.0e+03  -1.00e+00  2.121615263e+01   -1.011361719e+04  4.0e-02  0.13  
3   4.0e+01  9.7e-03  1.0e+03  -1.00e+00  9.052542078e+01   -9.961758548e+03  9.7e-03  0.13  
4   1.0e+01  2.4e-03  5.0e+02  -9.99e-01  3.666334519e+02   -9.360996350e+03  2.4e-03  0.14  
5   2.6e+00  6.3e-04  2.5e+02  -9.96e-01  1.411977590e+03   -7.110710012e+03  6.3e-04  0.14  
6   7.4e-01  1.8e-04  1.3e+02  -9.90e-01  4.880893881e+03   2.668326653e+02   1.8e-04  0.14  
7   2.9e-01  7.2e-05  8.4e+01  -9.80e-01  1.218738769e+04   1.556752270e+04   7.2e-05  0.14  
8   8.3e-02  2.0e-05  4.3e+01  -9.65e-01  4.156997801e+04   7.624242234e+04   2.0e-05  0.16  
9   3.3e-02  8.0e-06  2.6e+01  -9.19e-01  9.963655643e+04   1.921516644e+05   8.0e-06  0.16  
10  7.6e-03  1.8e-06  1.1e+01  -8.54e-01  3.462521309e+05   6.593094451e+05   1.8e-06  0.16  
11  1.6e-03  3.9e-07  2.5e+00  -4.81e-01  6.983636037e+05   1.101323345e+06   3.9e-07  0.16  
12  4.3e-04  1.1e-07  2.5e-01  8.46e-01   3.710592081e+05   4.260791013e+05   1.1e-07  0.17  
13  5.4e-05  1.3e-08  1.1e-02  9.97e-01   7.308564949e+04   7.991520525e+04   1.3e-08  0.17  
14  1.2e-07  3.0e-11  1.3e-06  1.00e+00   2.345063959e+02   2.523567123e+02   3.0e-11  0.17  
15  2.7e-09  6.6e-13  4.2e-09  1.00e+00   5.597032844e+00   5.981197787e+00   6.6e-13  0.17  
16  7.1e-10  1.7e-13  5.6e-10  1.00e+00   2.036489993e+00   2.139314120e+00   1.7e-13  0.19  
17  2.4e-10  5.8e-14  1.1e-10  1.01e+00   8.948401946e-01   9.297779685e-01   5.8e-14  0.19  
18  1.5e-10  3.5e-14  5.3e-11  1.01e+00   5.842215125e-01   6.059586407e-01   3.5e-14  0.20  
19  6.9e-11  1.7e-14  1.8e-11  1.00e+00   2.644252912e-01   2.750672154e-01   1.7e-14  0.20  
20  3.6e-11  8.8e-15  6.7e-12  1.01e+00   1.904933357e-01   1.960733133e-01   8.8e-15  0.22  
21  6.7e-11  3.4e-15  1.6e-12  1.01e+00   9.754699770e-02   9.971425118e-02   3.4e-15  0.22  
22  5.8e-11  3.4e-15  1.6e-12  9.94e-01   9.739141388e-02   9.955504312e-02   3.4e-15  0.22  
23  6.9e-11  3.4e-15  1.6e-12  9.93e-01   9.731044006e-02   9.947225024e-02   3.4e-15  0.23  
24  7.4e-11  3.4e-15  1.6e-12  9.99e-01   9.727358214e-02   9.943449754e-02   3.4e-15  0.23  
25  7.4e-11  3.4e-15  1.6e-12  9.95e-01   9.726366501e-02   9.942435650e-02   3.4e-15  0.25  
26  7.4e-11  3.4e-15  1.6e-12  9.94e-01   9.726366501e-02   9.942435650e-02   3.4e-15  0.25  
27  7.4e-11  3.4e-15  1.6e-12  9.94e-01   9.726366501e-02   9.942435650e-02   3.4e-15  0.27  
Optimizer terminated. Time: 0.33    


Interior-point solution summary
  Problem status  : UNKNOWN
  Solution status : UNKNOWN
  Primal.  obj: 9.7263665015e-02    nrm: 1e+07    Viol.  con: 6e-05    var: 0e+00    cones: 1e-09  
  Dual.    obj: 9.9424356596e-02    nrm: 6e-01    Viol.  con: 0e+00    var: 2e-12    cones: 0e+00  
Optimizer summary
  Optimizer                 -                        time: 0.33    
    Interior-point          - iterations : 28        time: 0.28    
      Basis identification  -                        time: 0.00    
        Primal              - iterations : 0         time: 0.00    
        Dual                - iterations : 0         time: 0.00    
        Clean primal        - iterations : 0         time: 0.00    
        Clean dual          - iterations : 0         time: 0.00    
    Simplex                 -                        time: 0.00    
      Primal simplex        - iterations : 0         time: 0.00    
      Dual simplex          - iterations : 0         time: 0.00    
    Mixed integer           - relaxations: 0         time: 0.00    

------------------------------------------------------------
Status: Inaccurate/Solved
Optimal value (cvx_optval): -0.0972637

That’s not a great solution, but then you get the warnings, and the norm of the primal is quite high. We don’t know how your problem looks numerically. Maybe you have a very wide range of coefficients or something else that should be improved. See https://docs.mosek.com/9.1/toolbox/debugging-numerical.html. If that doesn’t help contact Mosek support with a task file (cvx_solver_settings('write', 'dump.task.gz')).