Complementarity constraint

[blegat]

We discussed with @frapac how to create a complementarity constraint type.
One idea would be to add a type

struct Complementarity <: MOI.AbstractVectorFunction end
MOI.dimension(::Complementarity) = 2

And the constraint [x, y]-in-Complementarity would constraint x, y to be nonnegative and x * y to be zero.

[odow]

@xhub, @chkwon: you probably have some views on this.

At the JuMP level, it'd be nice to see something like this

@constraint(model, 0 <= F(x) ⟂ x >= 0)

[frapac]

The idea is to load directly complementarity constraints to the solver if it supports them.
For instance:
https://github.com/JuliaOpt/KNITRO.jl/blob/master/src/kn_constraints.jl#L362-L378

As far as I understand, Complementarity.jl reformulates complementarity constraints as NLP constraints. Once NLP would become first class in MOI, we may want to write bridges to reformulate complementarity constraints as NLP ones.

[mlubin]

SGTM. "How do I add complementarity constraints to JuMP?" was one of the original motivators for MOI.

[chkwon]

This will be great. In fact, I'm not so happy with the reformulations strategy currently in place with Complementarity.jl. A solver-native approach would be awesome.

By the way, I have used something like the following in Complementarity.jl

@complements(m, 0 <= F(x),   x >= 0)

A special symbol ⟂ is nice, but also we need a traditional ASCII approach as well. Perhaps:

@constraint(m, 0 <= F(x) complements   x >= 0)

[matbesancon]

it may be reformulated differently in MILP-based solvers, or solvers with better support for integrality than non-linear constraints, using SOS1 or indicator constraints

[xhub]

Sorry for giving my input so late, was too busy with ICCOPT 2019, and then I first got myself more familiar with MOI by porting my solver ReSHOP to MOI.

TLDR I would like to push for supporting more than just "complementarity constraint" (CC) like 0 <= F(x) ⟂ x >= 0, but LCP, MCP, VI (variational Inequalities), QVI (quasi-variational inequalities), MPCC and MPEC.

I started to type a design document 2 years ago for MBP, but right at the time the transition to MOI was started, so I will update it. I get the feeling that more problem classes can be handle through MOI.

I have a lightly different take on the kind of relations to support and the UI, which I will try to motivate by some perceived shortcomings of the proposed syntax. There is quite a few occurrences of the VI lingo, I hope it is not too off-putting.

Discussion on the "natural interface" @constraint(model, 0 <= F(x) ⟂ x >= 0)

With complementarity constraint, a salient point is correctness of the formulation, and the drawbacks with the kind of a "natural interface" like @constraint(model, 0 <= F(x) ⟂ x >= 0) include the following :

1. Nonsensical specification:@constraint(model, 0 == F(x) ⟂ x >= 0) is not a complementarity constraint.

2. Incorrect specification

@constraint(model, 0 <= F(x) ⟂ x >= 0)
@constraint(model, x <= 3)

The feasible set for x defines a rectangular VI (or MCP à la Ferris): F(x) ⟂ 0 <= x <= 3
Which case by case, means that:

• F(x) >= 0 whenever x = 0
• F(x) = 0 whenever 0 < x < 3
• F(x) <= 0 whenever x = 3

3. Inconsistent syntax: For a NCP, we would have a signature like @constraint(model, 0 <= F(x) ⟂ x >= 0) whereas that an MCP, it would be @constraint(m, F(x) ⟂ 0 <= x <= 3)

It's important to realize that in term of complementarity problems, the dual of set in which x belongs constraints is where F(x) belongs:
Hence, there is no need to specify twice the same information, also reducing the possible mistakes.

Here are 2 proposed ways to input complementarity relations

A) Popular CC

• Allow definition of rectangular VI via @constraint(m, F(x) ⟂ lb <= x <= ub), via ub and lb taken as +Inf, -Inf if not defined.
• An error is thrown if the user does something that breaks the type of constraint that he entered. Basically after the definition @constraint(m, 0 F(x) ⟂ lb <= x <= ub) sets x in stone, and no other constraints on x is allowed.
• If x was constrained before the CC was added, an error must be thrown.

Most of the work is in enforcing the "freezing" of x. Would it be done by JuMP, or left for the solver to be done? It being in JuMP would be nice to avoid duplicating this feature. On the MOI side, I am also not sure how to get from a VariableIndex all the ConstraintIndex where this variable is present. Maybe some solvers have this check internally.

B) General purpose complementarity specifications

• Provide an interface that support a much larger class of problem types: @complements(m, F, x)
  This would mean that the user doesn't specify directly a CC like A), but rather that there is a specific relation (named matching or complements) between a variable and a mapping.
• The actual type of constraint or problem, is not known until the optimize phase. The set where x belongs, is determined by looking at the the constraints involving only x. Then the proper problem class (NCP, MCP, VI, MPEC or MPCC) can be inferred. This is akin to an implicit splitting of the constraints in the model.

I favor B over A, but I have a limited view of the way other solvers (beyond ReSHOP and PATH) do support complementarity in models. I understand that both departs what most people are (unfortunately) taught about CC in the way that it is a relation where 0 <= F(x), x>=0, and F(x)^T x = 0, and that a majority of users may want to specify that kind of constraints. I hope that a good documentation/tutorial would clarify the point that there is no need for an explicit constraint on F.

MOI/JuMP rough implementation sketch for B

For B, at the MOI level, I would define a constraint [F, x]-in-generalized-equation, or [F, x]-in-VI that represent the generalized equation (another way of representing Variational Inequalities

where N_C is the normal cone from convex analysis,

Whenever C is closed convex cone, this is a (generalized) complementarity problem:

When C is the positive orthant or the whole space, one recovers "classical" complementarity condition.

At the JuMP level, I like the @complements(m, F, x) or @constraints(m, F ⟂ x), even @constraints(m, F ⟂ lb <= x <= ub) as a syntactic sugar that would translate into adding a constraint to x and setting the relationship between F and x. If the user specifies

@constraints(m, F ⟂ lb <= x)
@constraints(m, x <= ub)

This is understood as an MCP F ⟂ lb <= x <= ub

Quasi-Variational Inequalities (QVI)

For the Quasi-Variational Inequalities (QVI), one way to model them it to introduce a 3rd type of variable, y, such that the problem definition reads

At the solution, we require y = x. This would be covered at the MOI level by [F, x, y]-in-QVI. At the JuMP level, an interface of the kind @complements(F, x, y) or @QVI(F, x, y) could be considered

There are also opportunity for problem transformations (not an equivalence), MCP => NCP; VI => MCP; in addition to the aforementioned NCP=> binary/integer variables. In general, the first two should be avoided if the solver supports the original problem types, since the reformulated problem may not be processable. For instance, any VI over a compact set with a continuous functional has a solution (simple application of Browder FP theorem). The reformulation as a NCP may destroy that structure. This has been observed in practice over non-monotone problems.

And also relaxation MPCC => NCP, ..., have seen a lot of attentions in the past (NLPEC, work by S. Leyffer, ...).

In terms of solvers, AFAIK Baron does automatic detection. Knitro needs the mapping to be transformed into a defining equation for an auxiliary variable, which are then used in the CC definition.

[odow]

Thanks for the input @xhub!

I believe Michaels' preference is exactly your option B): complements(model, F, x), where F is a scalar (nonlinear) function and x is a variable, and that the bounds on x define the VI. (And for exactly the reasons you describe... (unsurprisingly).)

So we would define

struct Complements <: MOI.AbstractVectorSet end

Solvers (e.g., PATH.jl) would support

function MOI.add_constraint(model, ::MOI.VectorAffineFunction, ::MOI.Complements)

and throw an error is dimension(::MOI.VectorAffineFunction) != 2 or if the second dimension was not a straight variable.

Then, bridges would map people writing

@constraint(model, [y, x] in Complements())

into

@constraint(model, [1.0 y, x] in Complements())

JuMP can't (currently) handle

@constraint(model, [sin(y), x] in Complements())

because it doesn't have an AbstractScalarNonlinearFunction type. But this would be worked around by

@variable(model, z)
@NLconstraint(model, z == sin(y)
@constraint(model, [y, x] in Complements())

Handling quasi-variational inequalities seems harder, and a stretch goal. You could define

@constraint(model, [F(x), x, y] in QuasiVariationalInequality())

But the number of solvers supporting this seems minimal. It might be a good function/set combo for Path-VI to support until we work out the details, rather than add to MOI from the outset.

[xhub]

Would the following also work (inspired by SOS(1|2)) ?

struct Complements <: AbstractVectorSet
    variable::Union{MOI.VariableIndex,MOI.VectorOfVariables}
end

This would also allow for the use of vectors in the definition of CC.
I also like the semantics of having the variable and the set begin tied together.

The canonical signature would then be

function MOI.add_constraint(model, S, ::MOI.Complements) where S <: Union{MOI.VariableIndex, VOV, SF, VAF, VQF}

On the JuMP side, we would have

@constraint(model, y in Complements(x))
@constraint(model, x*x - 1 in Complements(x))

For the NLP workaround, it's a bit sad, I have read #846 . In EMP.jl I used some macro @NLvipair. Maybe we could have a nicer but hackisch @CCconstraint(m, func, x) which would translate into

@variable(m, z)
@NLconstraint(m, func == z)
@constraint(m, z in ComplementsNL(x))

where ComplementsNL <: Complements. I would push for a different signature for the VI/CP case.
Indeed, in that case if z were to stay as a variable, it would have to be matched with the constraint z = sin(x). For VI/CP, all variables must be defined has matching some function.

Furthermore, the variable z should be completely eliminated from the VI/CP since it most likely makes the feasible set nonconvex. Having a dedicated Set would be nice. I understand that one may not want to

For MPCC, the situation may be different, since for instance KNITRO requires the CC constraints to be between 2 variables. It seems that they would not need to differentiate.

Finally for QVI, how about a

struct QuasiVariationalInequalitySet <: AbstractVectorSet
    x::Union{MOI.VariableIndex,MOI.VectorOfVariables}
    y::Union{MOI.VariableIndex,MOI.VectorOfVariables}
end

That would allow for

@constraint(model, F(x) in QuasiVariationalInequality(x, y))

[odow]

I started an interface to PATH here: https://github.com/odow/PATH.jl.

(It's intended to replace https://github.com/chkwon/PATHSolver.jl, which uses an additional C wrapper instead of going directly into the libpath shared object, and has some funky ways of dealing with callbacks which are an artifact of old Julia code.)

My next step is to write the MOI interface. Talking with Michael, we will support VectorAffineFunction -in- Complements and VectorQuadraticFunction-in-Complements. At the JuMP level, you. will write

@constraint(model, [i=1:2], [i * x[1] + x[2], x[i]] in Complements())

You could also make F(x) a quadratic.

This is a little ungainly, so I might try to implement

@constraint(model, [i = 1:2], i * x[1] + x[2] ⟂ x[i])

Alternatively, we might take advantage of jump-dev/JuMP.jl#2051, and write

@constraint(model, [i=1:2], complements(i * x[1] + x[2], x[i]))

In any case, the JuMP syntax is purely cosmetic. The MOI syntax is the one required by the solvers.

Would the following also work (inspired by SOS(1|2)) ?

No. One limitation with the current MOI implementation is that we really don't want the sets to contain variables because this violates a lot of the assumptions made by bridges etc.

For the NLP workaround, it's a bit sad

I agree. In regard to nonlinear: I retract my previous work-around because you're right, it might destroy convexity/monotonicity.

We won't be targeting nonlinear equilibrium problems until MOI has a ScalarNonlinearFunction and VectorNonlinearFunction function. Then it should be seamless to write nonlinear F(x).

For MPCC, the situation may be different, since for instance KNITRO requires the CC constraints to be between 2 variables.

Does KNITRO allow

@constraint(model, z == x + y)
@constraint(model, [z, x] in Complements())

If so, it can just support VectorOfVariables-in-Complements and throw an error if the dimension is not 2.

[chkwon]

@odow A direct Julia interface to PATH is great! When PATHSolver.jl was written, installing PATH libraries was a bit of headache. Hope it will work out smoothly this time. By the way, they will not allow the name PATH.jl to be registered in METADATA.jl.

[odow]

installing PATH libraries was a bit of headache.

Yes, hopefully this will just require installing libpathXX, rather than compiling the Standalone_C interface. I believe we can use https://github.com/ampl/pathlib/tree/master/lib directly, although I have been developing for the latest Path 5.0.0. Michael was going to look into the API differences.

By the way, they will not allow the name PATH.jl to be registered in METADATA.jl.

Ah. Good to know. I guess PATH is a little generic.

[frapac]

@odow I confirm that Knitro would allow your formulation.

[xhub]

@frapac Based on the comments in https://github.com/JuliaOpt/KNITRO.jl/blob/master/examples/mpec1.jl I infered that it even requires that formulation. Am I missing something?

[xhub]

No. One limitation with the current MOI implementation is that we really don't want the sets to contain variables because this violates a lot of the assumptions made by bridges etc.

Ok that's unfortunate. Is it something that is planned on being worked on?

Sorry if I sound pushy, it would be great to be able to write

@constraint(model, A*x + b ⟂ x)

Could we have a

struct ComplementsPairing{F,V} <: AbstractFunction
    func::F
    variable::V
end

Then we can allow for scalar and vector use of Complements. The MOI signature would then be

function MOI.add_constraint(model, ::MOI.ComplementsPairing{F,V}, ::MOI.Complements)