Conic constraint in GAMS

GAMS
1

Hi experts!

I just ran into some problems for writing conic constraints in GAMS. I was trying to test the second-order cone programming (SOCP) duality by the following problem from the Stanford course:

min 2 * x1 + x2 + x3
s.t. x1 + x2 + x3 = 1
sqr(x1) =g= sqr(x2) + sqr(x3)

And the dual of this SOCP is

max y
s.t. y + s1 = 2
y + s2 = 1
y + s3 = 1
sqr(s1) =g= sqr(s2) + sqr(s3)

Seems simple right? Here is the code I wrote for the two programs, and I used Mosek as the solver.

free variable
         x1, x2, x3;

free variable
         y, obj, obj2;

equations
         eq1
         eq2
         eq3
         eq10;

eq1..            obj =e= 2 * x1 + x2 + x3;
eq2..            x1 + x2 + x3 =e= 1;
eq3..            sqr(x1) =g= sqr(x2) + sqr(x3);

option qcp = mosek;

model test /eq1, eq2, eq3/;

solve test using qcp minimizing obj;

free variable
         s1, s2, s3;

equations
         eq4
         eq5
         eq6
         eq7
         eq8;

eq4..            obj2 =e= y;
eq5..            y + s1 =e= 2;
eq6..            y + s2 =e= 1;
eq7..            y + s3 =e= 1;
eq8..            sqr(s1) =g= sqr(s2) + sqr(s3);

model test2 /eq4, eq5, eq6, eq7, eq8/;

solve test2 using qcp maximizing obj2;

If you try this code, Mosek will give you a piece of information saying that

Return code - 1293 [MSK_RES_ERR_CON_Q_NOT_PSD]: The quadratic constraint matrix is not PSD

Here PSD denotes positive semi-definite. But what? I believe the coefficient matrices in my eq3 and eq8 are PSD!

However, when I changed the notation of eq3 and eq8 to:

eq3..            x1 =c= x2 + x3;

and

eq8..            s1 =c= s2 + s3;

They now run correctly! And they both give the correct optimal solution of square root 2, which supports the strong duality.

Then I moved to Gurobi and Cplex, then found that the =c= operator is only valid in Mosek. But I still heard that Gurobi and Cplex CAN solve SOCP problems. So my question is, how can I code correctly by using =g= instead of =c= for a SOCP problem?

Thanks in advance! :smiley:

Gabriel

2

Hi,

with recent versions you should get a message like this from Mosek:

The constraint ‘eq8’(3) is not convex. Q should be negative semidefinite for a constraint with finite lower bound. (1294)

Setting non-negative lower bounds on x1 and s1 resolves the problem.

free variable
         x1, x2, x3;

x1.lo = 0;
free variable
         y, obj, obj2;

equations
         eq1
         eq2
         eq3
         eq10;

eq1..            obj =e= 2 * x1 + x2 + x3;
eq2..            x1 + x2 + x3 =e= 1;
eq3..            sqr(x1) =g= sqr(x2) + sqr(x3);

model test /eq1, eq2, eq3/;

option qcp = mosek;  solve test using qcp minimizing obj;
option qcp = cplexd; solve test using qcp minimizing obj;
option qcp = gurobi; solve test using qcp minimizing obj;



free variable
         s1, s2, s3;
s1.lo = 0
equations
         eq4
         eq5
         eq6
         eq7
         eq8;

eq4..            obj2 =e= y;
eq5..            y + s1 =e= 2;
eq6..            y + s2 =e= 1;
eq7..            y + s3 =e= 1;
eq8..            sqr(s1) =g= sqr(s2) + sqr(s3);

model test2 /eq4, eq5, eq6, eq7, eq8/;

option qcp = mosek;  solve test2 using qcp maximizing obj2;
option qcp = cplexd; solve test2 using qcp maximizing obj2;
option qcp = gurobi; solve test2 using qcp maximizing obj2;

I hope this helps!

Fred

3

Thanks Fred for your reply! Yes you are right. If we put a non-negative lower bound everything works fine, and that is how a convex cone works. It is just a little bit confusing that MOSEK still works and does not report any non-convexity if we do NOT set this non-negative lower bound, and it seems like it set the lower bound intrinsically. Gurobi and CPLEXD, however, need this lower bound to be set manually to let solvers work.

Gabriel

4

For example, see the following strong duality demo for a rotated conic programming:

variable x1, x2, x3, x4, obj;

equations
         objec, e1, e2, e3, e4;

option qcp = mosek;

objec..          obj =e= 2*x1 + x2 + x3 - 3*x4;
e1..             x1 =l= 7;
e2..             2*x1 + x3 - x2 + 0.5*x4 =e= 4;
e3..             x1 + x2 - 2*x4 =l= 15;
e4..             x1 + x2 =c= x3 + x4;

model test /all/;

solve test using qcp minimizing obj;

variable obj2, y1, y2, z1, z2, z3, z4;

positive variable y1, y3;

equations
         objec2, ee1, ee2, ee3, ee4, ee5;

objec2..         obj2 =e= - 7 * y1 + 4 * y2 - 15 * y3;
ee1..            2 + y1 - 2 * y2 + y3 - z1 =g= 0;
ee2..            1 + y2 + y3 - z2 =g= 0;
ee3..            1 - y2 - z3 =e= 0;
ee5..            -3 - 0.5*y2 - 2*y3 - z4 =e= 0;
ee4..            z1 + z2 =c= z3 + z4;

model test2 /objec2, ee1, ee2, ee3, ee4, ee5/;

solve test2 using qcp maximizing obj2;

We do not need to set the non-negative bound for x1, x2, and z1, z2, when MOSEK still works well with the =C= operator.

Gabriel

5

When using the =c= interface the cone is given explicitly including the non-negativity constraints. When you write as a quadratic program, MOSEK will try to recognize the cone, but if the bounds are not right than it refuses to do so.

-Michael