CONEPROG: Precision Issue (Exitflag -7)

3 views (last 30 days)
Holger I. Meinhardt on 14 Sep 2021
Commented: Aykut Bulut on 16 Sep 2021
These days I tried to port a symmetric positive semi-definite quadratic programming problem to a conic optimization problem. So-far I was able to transcribe the problem into a cone program, however it seems that still some options tuning is necessary, since all the iterations end up with an exitflag -7. Before I shall present the problem specification as a conic optimization problem and a snippet of the evaluation process, I shortly present a result of an example that is close to the solution. Moreover, it is worth noting that some examples finished up with the correct solution despite the consecutive appearance of the exitflag -7. In addition, I followed the piece of advice of the documentation to change the solver, but this had no effect. The same happened while changing the tolerance value from up to .
Due to the following returned outcome of an example, I suppose that I am on the right track, though the result below is not correct
x =
8.4966 11.5004 5.5087 23.0781 50.0829
The returned result of the evaluation is close to the expected solution below
pk_bv =
8.5000 11.5000 5.5000 23.0833 50.0833
which let me suppose that I encounter some numerical issues that I cannot localize. In the hope that a different angle may solve the puzzle, I present now the parameters of the problem. Moreover, we have , it should hold with equality.
The subsequent options setting has been chosen:
opts = optimoptions("coneprog","OptimalityTolerance",1e-12,"ConstraintTolerance",1e-12,"LinearSolver","auto")
And the cone problem is specified as:
Q=E'*E; % Matrix Q is symmetric positive semi-definite. Matrix E is made up by the numbers -1,0,1, it is not square.
b=-E'*a;
Aineq=[];
bineq=[];
c=[b(:);1]; % objective function
Asc = sqrtm(Q);
Asc((end+1),(end+1)) = 1;
d = [zeros(size(b(:)));1];
bsc = zeros(size(d));
gamma = -1;
socConstraints = secondordercone(Asc,bsc,d,gamma);
Aeq=[E(m,:),0];
lb=[vi(:);-inf];
ub=[ra(:);inf];
[u,fval,exitflag,output,lambda]=coneprog(c,socConstraints,Aineq,bineq,Aeq,a(m),lb,ub,opts)
In this respect, it might be useful to compare the above problem with the symmetric positive semi-definite quadratic programming problem:
Furthermore, I want now to present some snippets of the evaluation process to observe the options and output structure settings.
>> x=cp_kernel(bv)
opts =
coneprog options:
Set properties:
ConstraintTolerance: 1.0000e-12
LinearSolver: 'auto'
OptimalityTolerance: 1.0000e-12
Default properties:
Display: 'final'
MaxIterations: 200
MaxTime: Inf
Search direction is small and current iterate is not within the specified constraint and/or optimality tolerance.
u =
1.0e+04 *
0.0008
0.0012
0.0005
0.0023
0.0050
5.2886
fval =
-5.2785e+04
exitflag =
-7
output =
struct with fields:
iterations: 30
primalfeasibility: 1.0206e-11
dualfeasibility: 6.6721e-07
dualitygap: 9.6139e-12
algorithm: 'interior-point'
linearsolver: 'augmented'
message: 'Search direction is small and current iterate is not within the specified constraint and/or optimality tolerance.'
lambda =
struct with fields:
eqlin: 94.9844
ineqlin: [0x1 double]
soc: 9.9244e-09
lower: [6x1 double]
upper: [6x1 double]
ans =
8.4783 11.8906 5.2411 23.1990 49.8576
.
.
.
.
.
exitflag =
-7
output =
struct with fields:
iterations: 28
primalfeasibility: 1.0126e-11
dualfeasibility: 2.3589e-09
dualitygap: 1.0380e-11
algorithm: 'interior-point'
linearsolver: 'augmented'
message: 'Search direction is small and current iterate is not within the specified constraint and/or optimality tolerance.'
lambda =
struct with fields:
eqlin: -6.4447
ineqlin: [0x1 double]
soc: 2.7285e-08
lower: [6x1 double]
upper: [6x1 double]
Warning: Probably no kernel point found!
> In cp_kernel>computePrk (line 215)
In cp_kernel (line 102)
x =
8.4966 11.5004 5.5087 23.0781 50.0829
My question is now: which options must be tuned to get the correct result? Or has the cone problem specification a major flaw, which I have overlooked?
Thank you very much in advance for the kind support and any suggestion of improvements.
EDIT 15 Sep 2021
Thanks to the kind and excellent support of Aykut Bulut I could solve the above precision issue by reformulating the cone problem. For those who are looking for some inspiration for their own work, I present the reformulated cone problem by
Aineq=[];
bineq=[];
c=[zeros(size(b(:)));1]; % objective function
Asc = E;
Asc((end+1),(end+1)) = 1;
d = [zeros(size(b(:)));1];
bsc=a';
bsc(end+1)=0;
gamma = 0;
socConstraints = secondordercone(Asc,bsc,d,gamma);
Aeq=[E(m,:),0];
lb=[-inf(n+1,1)];
ub=[inf(n+1,1)];
[u,fval,exitflag,output,lambda]=coneprog(c,socConstraints,Aineq,bineq,Aeq,a(m),lb,ub,opts);

Aykut Bulut on 14 Sep 2021
Hi Holger,
The method you use for conversion is usefull when you do not have factors of Q. In your case, you can skip computing Q and using sqrtm() function to factor it. I think some precision might be lost in these operations and skipping them might improve the numerics.
Objective function is First observe that we can drop the constant term in the objective function. It can be added back once the optimal solution is computed. Second, instead of minimizing , we can minimize (note that we dropped the exponent 2 and \frac{1}{2} coefficient, these problems have the same set of optimal solutions and effectively same). Next, instead of minimizing , we define a new variable t and force it to be greater than equal to in the constraints, and minimize t. After this, objective function becomes t, and we have the following additional constraint. I think solving this modified problem (together with other constraints) will yield more precise solutions. This new problem will have a different optimal objective value. However, once you have optimal x for the new problem, you can compute the optimal objective function value of the original problem by using optimal x values computed.
Aykut Bulut on 16 Sep 2021
Thanks for your feedback. I will inform the team about the examples in the documentation.

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!