Incorrect solution?

[svigerske]

Issue created by migration from Trac.

Original creator: EsbenL

Original creation time: 2016-12-06 19:59:05

Assignee: ferreau

Version: 3.2.0

Hi. I am playing around with qpOASES and try to solve a QP that I have written myself, but I looks like I dont get a correct solution when the number of constraints exceeds the number of variables in the QP. In attached Images you can see the QP that I am trying to solve. Below you can see my qpOASES setup.

I am looking forward to hear from you. (Also I am sorry if this is not the right place to post this ticket, and then please let me know)

int main( )
{
    USING_NAMESPACE_QPOASES

    /* Setup data of first QP. */
    real_t H[2*2] = { 1.0, 0.0, 0.0, 1.0 };
    real_t g[2] = { 0.0, 0.0 };

    // c1
//    real_t A[1*2] = { 1.0, -1.0 };
//    real_t lbA[1] = { -1.0 };
//    QProblem example( 2,1 );      // => correct answer

    // c1, c2
//    real_t A[2*2] = { 1.0, -1.0,    0.5, 1.0};
//    real_t lbA[2] = { -1.0,2.0 };
//    QProblem example( 2,2 );      // => correct answer

    // c1, c2, c3
//    real_t A[3*2] = { 1.0, -1.0,    0.5, 1.0,   0.0, 1.0};
//    real_t lbA[3] = { -1.0, 2.0, 2.5 };
//    QProblem example( 2,3 );      // => incorrect answer

    // c1, c2, c4
//    real_t A[3*2] = { 1.0, -1.0,    0.5, 1.0,   3.0, -1.0};
//    real_t lbA[3] = { -1.0, 2.0, 3.0 };
//    QProblem example( 2,3 );      // => incorrect answer

    // c3, c4
    real_t A[2*2] = { 0.0, 1.0,   3.0, -1.0};
    real_t lbA[2] = { 2.5, 3.0 };
    QProblem example( 2,2 );        // => correct answer

    Options options;
    example.setOptions( options );
    int_t nWSR = 10;
    example.init( H,g,A,0,0,lbA,0, nWSR );

    /* Get and print solution of first QP. */
    real_t xOpt[2];
    real_t yOpt[2+1];
    example.getPrimalSolution( xOpt );
    example.getDualSolution( yOpt );
    printf( "\nxOpt = [ %e, %e ];  yOpt = [ %e, %e, %e ];  objVal = %e\n\n", 
            xOpt[0],xOpt[1],yOpt[0],yOpt[1],yOpt[2],example.getObjVal() );

    return 0;
}

[svigerske]

Attachment Graph.jpg by EsbenL created at 2016-12-06 19:59:57

QP

[svigerske]

Attachment Graph.2.jpg by EsbenL created at 2016-12-06 20:00:36

Graph showing QP

[svigerske]

Attachment QP.jpg by EsbenL created at 2016-12-06 20:02:23

QP

[svigerske]

Comment by ferreau created at 2016-12-07 07:59:32

Thanks for reporting. I am going to have a look at the issue soon.

[svigerske]

Comment by ferreau created at 2016-12-07 07:59:32

Changing status from new to assigned.

[svigerske]

Comment by ferreau created at 2016-12-07 07:59:32

Changing type from user support to test request.

[svigerske]

Comment by ferreau created at 2017-03-29 08:40:37

Apologies for taking such a long time to come back to you. I have tried all the five QP problems in your code snippet, but found all of them where solved correctly. In particular, I added to following lines at the end of the main function:

    SolutionAnalysis analyser;
    printf( "maxKktViolation: %e\n",analyser.getKktViolation(&example) );

    return 0;

to compute the maximum violation of the KKT optimality conditions. This value was 0 up to machine precision for all five QP problems. In case you did this already, note that there was a bug in the computation of this quantity (function getKktViolation in Utils.cpp) in case the QP features an identity Hessian matrix. That bug has been fixed in revision 233 in the development trunk and is going to become part of release 3.2.1.

Please let me know in case you have further comments on that issue. In particular, I could not fully relate the images to the supposedly incorrect solutions.

[svigerske]

Comment by ferreau created at 2017-04-03 08:49:00

Resolution: worksForMe

[svigerske]

Comment by EsbenL created at 2017-04-06 18:24:03

Thanks for your answer, I really appreciate it. Still I cant get the solver to give me correct/expected results in regards to xOpt. The solver does not give any errors, and therefore it looks like the QP is solved, but the optimal solution is not what I expect. Let me walk you through my example in more details.

I will again use the QP described earlier .. so please have a look at the image “Graph.jpg” that visually shows the QP together with the QP-math in the other attached image “QP.jpg”. As you can see there are 4 constraints numbered 1,2,3 and 4. Now when I setup QPs with different combinations of these contraints I do not see results from the solver as i would have expected, when looking at the image. Below I will show you first a QP with a result that I agree is optimal, and hereafter a QP with a result that I do not agree is optimal.

QP with contraints 3 and 4:

    real_t A[2*2] = { 0.0, 1.0,   3.0, -1.0};
    real_t lbA[2] = { 2.5, 3.0 };
    QProblem example( 2,2 );

gives the following result: 
####################   qpOASES  --  QP NO.   1   #####################

    Iter   |    StepLength    |       Info       |   nFX   |   nAC    
 ----------+------------------+------------------+---------+--------- 
       0   |   2.255636e-07   |   ADD CON    1   |     1   |     1   
       1   |   1.250851e-07   |   ADD CON    0   |     0   |     2   
       2   |   1.000000e+00   |    QP SOLVED     |     0   |     2   

xOpt = [ 1.833333e+00, 2.500000e+00 ];  yOpt = [ 0.000000e+00, 0.000000e+00, 3.111111e+00 ];  objVal = 4.805556e+00

And when looking at the Image, I agree that this is indeed the optimal solution ( please have a look at the image with focus on xOpt just found ).

Now lets try to add one more constraint which gives us : QP with contraints 2, 3 and 4:

    real_t A[3*2] = { 0.5, 1.0,    0.0, 1.0,    3.0, -1.0};
    real_t lbA[3] = { 2.0, 2.5, 3.0 };
    QProblem example( 2,3 );

gives the following result: 
####################   qpOASES  --  QP NO.   1   #####################

    Iter   |    StepLength    |       Info       |   nFX   |   nAC    
 ----------+------------------+------------------+---------+--------- 
       0   |   2.311108e-07   |   ADD CON    2   |     1   |     1   
       1   |   9.251744e-09   |   ADD CON    0   |     0   |     2   
       2   |   1.684781e-01   |   ADD CON    1   |     0   |     2   
       3   |   1.000000e+00   |    QP SOLVED     |     0   |     2   

xOpt = [ 6.111111e-01, 2.500000e+00 ];  yOpt = [ 0.000000e+00, 0.000000e+00, 0.000000e+00 ];  objVal = 4.805556e+00

And when looking at the Image again, I do not agree that its the optimal point. As you can see from the Image, the newly added constraint (nr 2) will not affect the optimal solution and therefore I would expect the result to be the same as before. I do see that objVal found in both scenarios are the same, but I dont understand why xOpt are not the same.

This Is what I do not understand. Please let my know what you think.

Best Esben

[svigerske]

Comment by ferreau created at 2017-04-10 09:48:07

Thanks for walking me through your example. I agree that the second solution you posted in not correct. However, I cannot reproduce your results. Both C++ core and Matlab interface also give the correct solution in the second case (and the additional constraint indeed does not effect it). As the first entry of your xOpt matches an entry in the dual solution: could it be that you have not increased the size of array yOpt for your second example and you see the artefact of a memory access violation?

[svigerske]

Comment by EsbenL created at 2017-04-10 11:35:10

Thanks for your reply.

My example implementation looks like this:

int main( )
{
    USING_NAMESPACE_QPOASES

    /* Setup data of first QP. */
    real_t H[2*2] = { 1.0, 0.0, 0.0, 1.0 };
    real_t g[2] = { 0.0, 0.0 };

    //    // c2, c3, c4
    real_t A[3*2] = { 0.5, 1.0,    0.0, 1.0,    3.0, -1.0};
    real_t lbA[3] = { 2.0, 2.5, 3.0 };
    QProblem example( 2,3 );      // => correct solution : xOpt = [ -1.110223e-16, 2.500000e+00 ];  yOpt = [ 0.000000e+00, 0.000000e+00, 0.000000e+00 ];  objVal = 3.125000e+00

    Options options;
    options.printLevel = PL_MEDIUM;

    example.printProperties();

    example.setOptions( options );

    int_t nWSR = 10;
    example.init( H,g,A,NULL,NULL,lbA,NULL, nWSR );

    /* Get and print solution of first QP. */
    real_t xOpt[2];
    real_t yOpt[3];
    example.getPrimalSolution( xOpt );
    example.getDualSolution( yOpt );
    printf( "\nxOpt = [ %e, %e ];  yOpt = [ %e, %e, %e ];  objVal = %e\n\n", 
            xOpt[0],xOpt[1],yOpt[0],yOpt[1],yOpt[2],example.getObjVal() );

    SolutionAnalysis analyser;
    printf( "maxKktViolation: %e\n",analyser.getKktViolation(&example) );

    return 0;
}

Does it look correct to you ?

also when I compile it, I do get a few warnings :

Creating _myTest.o
clang: warning: optimization flag '-finline-functions' is not supported
In file included from _myTest.cpp:36:
In file included from /Users/esbenlundsager/Documents/StudioLundsager/Development/Optimization/qpOASES/qpOASES-3.2.0/include/qpOASES.hpp:62:
/Users/esbenlundsager/Documents/StudioLundsager/Development/Optimization/qpOASES/qpOASES-3.2.0/include/qpOASES/SQProblemSchur.hpp:196:23: warning: 
      'qpOASES::SQProblemSchur::setupAuxiliaryQP' hides overloaded virtual
      function [-Woverloaded-virtual]
                virtual returnValue setupAuxiliaryQP(   SymmetricMatrix ...
                                    ^
/Users/esbenlundsager/Documents/StudioLundsager/Development/Optimization/qpOASES/qpOASES-3.2.0/include/qpOASES/QProblem.hpp:848:23: note: 
      hidden overloaded virtual function 'qpOASES::QProblem::setupAuxiliaryQP'
      declared here: different number of parameters (2 vs 6)
                virtual returnValue setupAuxiliaryQP(   const Bounds* co...
                                    ^
1 warning generated.
Creating /Users/esbenlundsager/Documents/StudioLundsager/Development/Optimization/qpOASES/qpOASES-3.2.0/bin/_myTest
rm _myTest.o

But I am not sure what this means. you? Looking forward to hearing from you again. Best Esben

[svigerske]

Comment by ferreau created at 2017-04-10 11:43:21

It should read real_t yOpt[5]; as there are also two Langrange multiplier for the simple bounds (even if you did not specify them).

Regarding the warnings: compiler flag -finline-functions has been removed in the latest release 3.2.1. I am going to have a look at the other ones, which, howewer, is surely unrelated to the optimal solution issue.

[svigerske]

Comment by EsbenL created at 2017-04-10 12:29:36

wow then it works! thank you so much!