'drake gurobi cannot handle quadratic constraint

I am currently working to implement a simple quadratic constraint on an optimization problem. Gurobi's website says that quadratic constraints can be implemented. Does drake not have an interface with this constraint using Gurobi?

The code is as follows.

#include <Eigen/Dense>
#include <math.h>
#include <iostream>

#include "drake/solvers/mathematical_program.h"
#include "drake/solvers/solve.h"
#include "drake/solvers/gurobi_solver.h"
#include "drake/solvers/scs_solver.h"

using namespace std;
using namespace Eigen;
using namespace drake;

int main(){
    solvers::MathematicalProgram prog; 

    // Instantiate the decision variables
    solvers::VectorXDecisionVariable x = prog.NewContinuousVariables(2, "x");

    // Define constraints 
    for(int i = 0; i < 2; i++){
      prog.AddConstraint(x[i]*x[i] <= 2); //Replacing this with a linear constraint 
                                          //such as prog.AddLinearConstraint(-5 <= x[i] && x[i] <= 5); 
                                          //will work
    }

    // Define the cost
    MatrixXd Q(2,2);
    Q << 2, 1,
         1, 4;
    VectorXd b(2);
    b << 0, 3;
    double c = 1.0;
    prog.AddQuadraticCost(Q, b, c, x);

    solvers::GurobiSolver solver;
    cout << "Gurobi available? " << solver.is_enabled() << endl;

    auto result = solver.Solve(prog);
    cout << "Is optimization successful?" << result.is_success() << endl;
    cout << "Optimal x: " << result.GetSolution().transpose() << endl;
    cout << "solver is: " << result.get_solver_id().name() << endl;
    cout << "computation time is: " << result.get_solver_details<solvers::GurobiSolver>().optimizer_time;

    return 0;
}
drake


Solution 1:[1]

Drake currently (intentionally) doesn't support the quadratic constraint. For better solver performance, I would recommend to reformulate the optimization problem without quadratic constraint, but with second-order cone constraints.

First for a quadratic constraint

x?x <= a²

This can be regarded as a Lorentz cone constraint

(y, x) is in the Lorentz cone
y = a

In Drake, you could add this constraint as

prog.AddLorentzConeConstraint((drake::VectorX<drake::symbolic::Expression>(x.rows() + 1) << a, x).finished());

For the more general quadratic constraint 0.5x?Qx + b?x + c<=0, we also provide AddQuadraticAsRotatedLorentzConeConstraint which imposes a rotated Lorentz cone constraint.

prog.AddQuadraticAsRotatedLorentzConeConstraint(Q, b, c, x);

For the quadratic cost

min 0.5x?Qx + b?x + c

you could also reformulate it as a linear cost with a rotated Lorentz cone constraint

min z
z >= 0.5x?Qx + b?x + c

In Drake, the code is

z = prog.NewContinuousVariables<1>()(0);
prog.AddLinearCost(z);
prog.AddRotatedLorentzConeConstraint(symbolic::Expression(z), symbolic::Expression(1), 0.5 * x.dot(Q * x) + b.dot(x) + c);

The reason why we would prefer Lorentz cone constraint, instead of a quadratic-constrained-quadratic-program (QCQP), is that with Lorentz cone constraint, you end up with a second-order conic optimization problem (SOCP), this optimization problem has a lot of nice features, in terms of its dual form. You could read the paper Alizadeh's paper for more details. And Mosek also recommends SOCP over QCQP https://docs.mosek.com/9.2/rmosek/prob-def-quadratic.html#a-recommendation

Solution 2:[2]

Possibly the (as-of-yet unimplemented) feature request at https://github.com/RobotLocomotion/drake/issues/6341 is what you're looking for?

Sources

This article follows the attribution requirements of Stack Overflow and is licensed under CC BY-SA 3.0.

Source: Stack Overflow

Solution Source
Solution 1
Solution 2 jwnimmer-tri

Related Questions