Add PSD constraint using Burer-Monteiro formulation

[hongkai-dai]

We can formulate the convex constraint X is PSD with a non-convex constraint using the Burer-Monteiro formulation

X = YYᵀ

This can be useful when we have non-convex optimization problem (for example, a bilinear matrix inequality), and we can solve the problem through gradient based nonlinear optimization (like SNOPT).

To support this feature, I propose that we add a class

/**
Takes X (the lower diagonal part) and Y (the lower diagonal part), and evaluates X - YYᵀ (the lower diagonal part).
*/
class BurerMonteiroPsdConstraint: public Constraint {
};

and the function

std::pair<Binding<BurerMonteiroPsdConstraint>, VectorX<symbolic::Variable>> MathematicalProgram::AddBurerMonteiroPsd(const MatrixX<symbolic::Variable>& X);

which also adds the slack variable Y into MathematicalProgram.

WDYT @AlexandreAmice @RussTedrake

[AlexandreAmice]

It would be nice to be able to support Burer-Monteiro, though I have no idea how SNOPT would do on that problem.

I am a little weary of adding another class though as it seems unnecessary. Would something like the design pattern for LorentzCone where we have a PsdConstraint::EvalType::kBurerMonteiro be better? That way there is only one concept for PsdConstraint. The main challenge here is how to allow the user to specify the rank of Y.

This way, if the user wants to specify a non-linear formulation and use a non-linear solver they can, but if we were to say support SDPLR or HALLaR, we don't need to explicitly handle both the Burer-MonteiroConstraint and PsdConstraint separately, there would only be one entry point.

[hongkai-dai]

I had the same concern of creating BurerMonteiroPsdConstraint. But adding a flag PsdConstraint::EvalType:kBurerMonteiro also has its pitfalls.

Evaluating the Burer-Monteiro formulation requires knowing not only X, but also the slack Y variable. This is different from the current PositiveSemidefiniteConstraint, which only requires knowing X. Hence depending on the EvalType, the number of variables for this constraint will change.

Also if we call AddPositiveSemidefiniteConstraint(psd_constraint, vars), what should be included this vars if EvalType == kBurerMonteiro? I think to keep the signature consistent with other EvalType, vars should only include X. And this function would add the slack Y internally.

[AlexandreAmice]

I think what you would do here is have Y be an internal member of PsdConstraint. It would NOT be added to the mathematical program, it would only be passed onto the solver. Effectively, there is no difference between X and Y. Each can be recovered by the other and so we should not keep two copies around.

In fact, I think if AddPositiveSemidefiniteConstraint(psd_constraint, vars, EvalType:kBurerMonteiro) then vars should be interpretted as the Y variable and only the Y variable should be stored to avoid the n*n storage.