Set inequality syntax suggestion

[odow]

See #3769 (comment)

Problem

Consider:

model = Model();
@variable(model, x[1:2, 1:2])
@constraint(model, x >= 0)

Does this mean:

@constraint(model, x >= 0, Nonnegatives())

or

@constraint(model, x >= 0, PSDCone())

Because of this ambiguity, we currently disallow the x >= 0 syntax for matrix sets.

Suggestion

@LebedevRI suggests that we introduce the ⪰ operator with the following behavior:

@constraint(model, x ⪰ y, set)
# lowers to
@constraint(model, x - y in set)

Examples

@constraint(model, x >= 0)
# lowers to
@constraint(model, x ⪰ 0, Nonnegatives())

@constraint(model, x <= 0)
# lowers to
@constraint(model, x ⪰ 0, Nonpositives())

@constraint(model, x ⪰ 0)   
# ERROR: missing a set that defines the partial order

@constraint(model, x >= 0, set)  
# lowers to
@constraint(model, x ⪰ 0, set)

Tradeoffs

The biggest pro is that it clarifies the distinction between ≥ and ⪰.

This is also the biggest con: I worry that users will not make the distinction, and people will write x >= 0 expecting that to be a PSD constraint when it is really a nonnegative constraint.

Currenntly, they need to opt-in by passing a set to avoid the ambiguity.

[LebedevRI]

@odow thank you for posting this for disscussion!

My main concern is syntax divergence. With the current behavior,
given current JuliaLang/julia#3766 + JuliaLang/julia#3769, while @constraint(model, x >= y) will work for scalars/vectors,
but if you've got matrices, you will need to use a completely different syntax:
@constraint(model, x >= y, Nonnegatives) or @constraint(model, x - y, Nonnegatives),
which requires actual reasoning about which of the sets should be used there,
wherein with @constraint(model, x >= y) syntax it just works.

I agree that this is very much not a clear-cut choice,
though it would seem to me that making such a change
would allow for a much more natural-looking models.

[odow]

My answer to:

while @constraint(model, x >= y) will work for scalars/vectors, but if you've got matrices, you will need to use a completely different syntax

is "yes, because of the ambiguity that users will experience, particularly those coming from YALMIP, where >= defaults to PSDCone: https://jump.dev/JuMP.jl/dev/tutorials/transitioning/transitioning_from_matlab"

There is no ambiguity for scalars and vectors.

[LebedevRI]

My answer to:

while @constraint(model, x >= y) will work for scalars/vectors, but if you've got matrices, you will need to use a completely different syntax

is "yes, because of the ambiguity that users will experience, particularly those coming from YALMIP, where >= defaults to PSDCone: https://jump.dev/JuMP.jl/dev/tutorials/transitioning/transitioning_from_matlab"

Ah, i did not know that. But then, is it not reasonable to expect them to read the documentation?
They'll also be able to notice the problem by actually looking at the built constraint (by printing it/model).
It seems unfortunate that said short transitonary/learning period expirience
should make better experience for everyone else be unachievable.

[blegat]

I think x >= y between scalars is clear, x >= y between vectors is slightly ambiguous because it does not say what is the order between vectors but I don't think there is too much risk that someone would misinterpret this and think that another cone than the positive orthant would be used by default.
For matrix cones, it's an inequality with no cone specified could be interpreted differently by different people so I don't think it would be helpful. Even for users that have read the docs, using this syntax would make their model have the potential of being misinterpreted.
For many people, >= is the only way to compare things in ascii so adding subtle difference with unicode character is going to be a source of confusion.

[LebedevRI]

Copying more of #3769 (comment) here, for clarity:
The proposed behavior is:

@constraint(model, x \succeq y)
# error: specify one of the conic sets

@constraint(model, x >= y, <... Some conic set ...>)
# warning: use \succeq instead

@constraint(model, x \succeq y, <... Some non-conic set ...>)
# warning: use >= instead

@constraint(model, x >= y)
# becomes @constraint(model, x - y, Nonnegatives)

@constraint(model, x <= y)
# becomes @constraint(model, x - y, Nonpositives)

I.e. if you use >=/<=, and don't specify a set, or specify a non-conic set, you get vector behavior,
if you use >=/<=, and specify a conic set, you get a warning to use \succeq instead,
and if you use \succeq, you must specify a set, and if it's not a conic one, you get a warning.

So really, the only thing one needs to remember is to use \succeq when meaning conic constraint.
Is this really too much to ask of users?

[odow]

1. We could do this
2. We can't do this. That is breaking. We don't want to start warning for commonly used existing code.
3. We can't do this. There is no "conic"/"non-conic" distinction in MOI
4. We could do this, but there is the ambiguity question
5. We could do this, but there is the ambiguity question

Is this really too much to ask of users?

Experience has told me: yes. Especially for novice users who are most at risk of making this mistake. (Not to pick on one example, but here is a recent one that I remember: https://discourse.julialang.org/t/lmi-optimization-problem/112218. They ended up trying .<= 0 instead of <= 0, PSDCone() just to get something that didn't error.)

I'm okay adding the \succeq operator as a new opt-in feature. I don't want to change existing behavior.

[blegat]

What do you mean by @constraint(model, x \succeq y, <... Some non-conic set ...>) ? What is an example of non-conic set ? Do you mean Nonnegatives is not a cone ? The issue is you can also see Nonnegatives as just one cone among the others and it's not really the first one you would think of when thinking of matrices.

[jd-foster]

This discussion is reminiscent of JuliaLang/julia#44358 where a breaking change in Base Julia was made to avoid element vs. vector division confusion.

[LebedevRI]

What do you mean by @constraint(model, x \succeq y, <... Some non-conic set ...>) ?

I specifically mean PSDCone/HermmitianPSDCone, yes, which are, per previous comments,
what newcomers may think what the matrix x >= y may mean.

What is an example of non-conic set ?

I only mean Zeros/Nonnegatives/Nonpositives by that.

This discussion is reminiscent of JuliaLang/julia#44358 where a breaking change in Base Julia was made to avoid element vs. vector division confusion.

Right, though here the behavior would be consistent (and in fact,
remove an inconsistency that the matrix and non-matrix arguments
require different syntax), and well-documented. The fear really is that
newcomers may not spare time to actually read what the docs say
and write one thing thinking it does something completely different...

[odow]

Thinking on this over the last few days, I don't even particularly want to add \succeq as an opt-in syntactic sugar for what we already have. It opens up the question of "how is it different from >=" and our answer will be that it isn't.

What we currently have is Pareto optimal on a set of tradeoffs. We could make it better for some people at the expense of others, and it isn't obvious that moving one particular way is preferable.

[LebedevRI]

@blegat

For matrix cones, it's an inequality with no cone specified could be interpreted differently by different people so I don't think it would be helpful. Even for users that have read the docs, using this syntax would make their model have the potential of being misinterpreted.

1. That's why i'm suggesting \succeq as the "only" method to introduce conic constraints
   (as per user-facing documentation, i get the part about not breaking existing models).
   With that disambiguation >= becomes available to be used for the $subj syntax sugar.
2. Isn't that confusion essentially already there? Given just x >= y, you can't know if they are matrices or not?
   I don't get this argument.

@odow

Thinking on this over the last few days, I don't even particularly want to add \succeq as an opt-in syntactic sugar for what we already have. It opens up the question of "how is it different from >=" and our answer will be that it isn't.

Oh come on...

">= syntax is used to perform element-wise comparisons, regardless the type of the variable,
as long as the shapes match, while the \succeq syntax is used when
a) you are operating on a matrix, and b) you are adding a explicitly specified conic constraint.
While >= syntax can also be used to add such constraints, that is a legacy interface
that is retained only as a stop-gap measure to support legacy models,
it may be dropped in future major version."

What we currently have is Pareto optimal on a set of tradeoffs. We could make it better for some people at the expense of others, and it isn't obvious that moving one particular way is preferable.

Perhaps i'm missing some context that was communicated somewhere elsewhere,
but i'm still a bit confused what those exact groups are.

From where i'm sitting, i wouldn't expect those who have previously encountered
this issue (using/thinking >= 0 to mean `>= 0, PSDCone()``) before, to completely forget it the next day.
Is that an unreasonable expectation? That's the delimiter i thought is in play,
and then i would say one group is disproportionally larger than the other one.

Experience has told me: yes. Especially for novice users who are most at risk of making this mistake. (Not to pick on one example, but here is a recent one that I remember: https://discourse.julialang.org/t/lmi-optimization-problem/112218. They ended up trying .<= 0 instead of <= 0, PSDCone() just to get something that didn't error.)

As someone who had to deal with users before, to me that tells that even with the current behavior
(a diagnostic "did you mean to use PSDCone()?"), there will always be cases where a completely wrong thing
will be done, falsely hoping it does something else. Trying to avoid that at all cost seems impractical, IMO.
It will always happen regardless.

[blegat]

The current design is pareto optimal with respect to the readability of models. If you see x >= y alone and you know they are vectors (ok if you don't know their shape it could be ambigious, that's why we didn't directly add it). If you see A >= B alone and A and B are symmetric square matrices, some people will see this a matrices, because it is a much more natural order between symmetric square matrices. No amount of documentation can solve that. Some people may loose hours debugging before finally figuring out that this actually means elementwise inequality.

The fear really is that newcomers may not spare time to actually read what the docs say and write one thing thinking it does something completely different...

Yes, or use JuMP again after some time so don't read the doc again, or know the distinction but made a mistake by using the wrong unicode character and don't notice it because you need to zoom in to see the difference.

With the cone specified (including Nonnegatives as a cone here), it's always explicit and clear. You need more keystrokes but also increase the readability of your model so each of these keystrokes is a win IMO. JuMP tries to be a modeling language that works for Constraint Programming, Conic Programming, Linear Programming, Non-Linear Programming, ... so we can't use a syntax that will feel natural for LP but that will make people from other fields struggle

[LebedevRI]

Okay, let's backtrack.

My syntax complaint is: it seems like having to manually go between the x >= y / x <= y syntax, and
@constraint(model, x - y in Nonnegatives) / @constraint(model, x - y in Nonpositives) is lacking sugar.
While, clearly, the suggestion of allowing @constraint(model, x >= y) to just work for matrices
is not receiving a warm welcome, there is at least one alternative:

Introduce a 'phony' JuMP set, ElementWiseComparison (name TBD),
which is immediately lowered into either Nonnegatives or Nonpositives, depending on the comparison:

@constraint(model, x >= y, ElementWiseComparison)
   =>
@constraint(model, x - y in Nonnegatives)

@constraint(model, x <= y, ElementWiseComparison)
   =>
@constraint(model, x - y in Nonpositives)

@constraint(model, x >= y)
# error: <...>

@constraint(model, x <= y)
# error: <...>

@constraint(model, x <= y, PSDCone)
# error: use `>=`

Hopefully surely that is a reasonable acceptable middle ground solution? :)

[odow]

First, note that if you want element wise inequality, you have always been able to write:

@constraint(model, x .<= y)

This is perhaps the simplest and most explicit way of writing it down. The only downside is that in the symmetric matrix case, it does not exploit symmetry. (But perhaps we could make it.)

Nonnegatives is pretty much your ElementWiseComparison, except that we would have:

@constraint(model, x <= y, Nonnegatives())
   =>
@constraint(model, y - x in Nonnegatives())

Because of the slight mis-match, we currently error:
https://github.com/jump-dev/JuMP.jl/pull/3766/files#diff-8cd8322515b4f61e46a3296d7ca78445d1eb335270416ee1dec6dd831b5e8ad7R686-R695

We can't do this because it breaks existing syntax:

@constraint(model, x <= y, PSDCone)
# error: use `>=`

[LebedevRI]

First, note that if you want element wise inequality, you have always been able to write:

@constraint(model, x .<= y)

I'm very aware of that, yes.

This is perhaps the simplest and most explicit way of writing it down. The only downside is that in the symmetric matrix case, it does not exploit symmetry. (But perhaps we could make it.)

Yes, that might be good. It seems not optimal to create duplicate constraints,
but then such deduplication requires mode code.

Nonnegatives is pretty much your ElementWiseComparison, except that we would have:

@constraint(model, x <= y, Nonnegatives())
   =>
@constraint(model, y - x in Nonnegatives())

Because of the slight mis-match, we currently error: https://github.com/jump-dev/JuMP.jl/pull/3766/files#diff-8cd8322515b4f61e46a3296d7ca78445d1eb335270416ee1dec6dd831b5e8ad7R686-R695

We can't do this because it breaks existing syntax:

@constraint(model, x <= y, PSDCone)
# error: use `>=`

I'm getting increasingly more and more confused.
Is it, at least, clear as to why i'm asking for ElementWiseComparison syntax sugar(1)?
There is either some very basic misunderstanding going on, or something else less friendly.

We can't do this because it breaks existing syntax:

I'm afraid i don't quite follow. Firstly, it does not break existing syntax,
because the lines you are referencing are in an open unmerged pull request.

julia> @constraint(model, x >= y, PSDCone())
[x[1,1] - y[1,1]  x[1,2] - y[1,2];
 x[2,1] - y[2,1]  x[2,2] - y[2,2]] ∈ PSDCone()

julia> @constraint(model, x <= y, PSDCone())
[-x[1,1] + y[1,1]  -x[1,2] + y[1,2];
 -x[2,1] + y[2,1]  -x[2,2] + y[2,2]] ∈ PSDCone()

secondly, i'm specifically asking about ElementWiseComparison,
so i'm not following what PSDCone has to do with it.

It doesn't break syntax for ElementWiseComparison because
we don't syntax for ElementWiseComparison because
we don't have ElementWiseComparison in the first place.

Again, i don't actually care how exactly ElementWiseComparison is expanded,
only that @constraint(model, x_matrix >= y_matrix, ElementWiseComparison)
does the same thing as @constraint(model, x_matrix .>= y_matrix),
and that @constraint(model, x_matrix <= y_matrix, ElementWiseComparison)
does the same thing as @constraint(model, x_matrix .<= y_matrix)
(modulo vectorization/symmetry)...

Perhaps i'm on the wrong end of the bell curve and missing a few decades of LP expirience,
can you please help me understand why, specifically, that can not be done?

[odow]

Let's make this discussion more concrete.

What we currently support

Here is what we currently support.

Vectors

julia> model = Model();

julia> @variable(model, x[1:2])
2-element Vector{VariableRef}:
 x[1]
 x[2]

julia> y = [1, 2]
2-element Vector{Int64}:
 1
 2

julia> @constraint(model, x .>= y)
2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.GreaterThan{Float64}}, ScalarShape}}:
 x[1] ≥ 1
 x[2] ≥ 2

julia> @constraint(model, x .<= y)
2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}, ScalarShape}}:
 x[1] ≤ 1
 x[2] ≤ 2

julia> @constraint(model, x >= y)
[x[1] - 1, x[2] - 2] ∈ MathOptInterface.Nonnegatives(2)

julia> @constraint(model, x <= y)
[x[1] - 1, x[2] - 2] ∈ MathOptInterface.Nonpositives(2)

julia> @constraint(model, x >= y, SecondOrderCone())
[x[1] - 1, x[2] - 2] ∈ MathOptInterface.SecondOrderCone(2)

julia> @constraint(model, x <= y, SecondOrderCone())
[-x[1] + 1, -x[2] + 2] ∈ MathOptInterface.SecondOrderCone(2)

We don't currently support Nonnegatives et al:

julia> @constraint(model, x <= y, Nonnegatives())
ERROR: At REPL[21]:1: `@constraint(model, x <= y, Nonnegatives())`: Unrecognized constraint building format. Tried to invoke `build_constraint(error, AffExpr[x[1] - 1, x[2] - 2], MathOptInterface.LessThan{Bool}(false), Nonnegatives())`, but no such method exists. This is due to specifying an unrecognized function, constraint set, and/or extra positional/keyword arguments.
[...]

julia> @constraint(model, x >= y, Nonnegatives())
ERROR: At REPL[22]:1: `@constraint(model, x >= y, Nonnegatives())`: Unrecognized constraint building format. Tried to invoke `build_constraint(error, AffExpr[x[1] - 1, x[2] - 2], MathOptInterface.GreaterThan{Bool}(false), Nonnegatives())`, but no such method exists. This is due to specifying an unrecognized function, constraint set, and/or extra positional/keyword arguments.
[...]

Matrices

julia> model = Model();

julia> @variable(model, x[1:2, 1:2])
2×2 Matrix{VariableRef}:
 x[1,1]  x[1,2]
 x[2,1]  x[2,2]

julia> y = [1 2; 2 3]
2×2 Matrix{Int64}:
 1  2
 2  3

julia> @constraint(model, x .>= y)
2×2 Matrix{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.GreaterThan{Float64}}, ScalarShape}}:
 x[1,1] ≥ 1  x[1,2] ≥ 2
 x[2,1] ≥ 2  x[2,2] ≥ 3

julia> @constraint(model, x .<= y)
2×2 Matrix{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}, ScalarShape}}:
 x[1,1] ≤ 1  x[1,2] ≤ 2
 x[2,1] ≤ 2  x[2,2] ≤ 3

julia> @constraint(model, x >= y, PSDCone())
[x[1,1] - 1  x[1,2] - 2;
 x[2,1] - 2  x[2,2] - 3] ∈ PSDCone()

julia> @constraint(model, x <= y, PSDCone())
[-x[1,1] + 1  -x[1,2] + 2;
 -x[2,1] + 2  -x[2,2] + 3] ∈ PSDCone()

We don't support no cone or Nonnegatives:

julia> @constraint(model, x >= y)
ERROR: At REPL[32]:1: `@constraint(model, x >= y)`: Unsupported matrix in vector-valued set. Did you mean to use the broadcasting syntax `.>=` instead? Alternatively, perhaps you are missing a set argument like `@constraint(model, X >= 0, PSDCone())` or `@constraint(model, X >= 0, HermmitianPSDCone())`.
Stacktrace:
[...]

julia> @constraint(model, x <= y)
ERROR: At REPL[33]:1: `@constraint(model, x <= y)`: Unsupported matrix in vector-valued set. Did you mean to use the broadcasting syntax `.<=` instead? Alternatively, perhaps you are missing a set argument like `@constraint(model, X <= 0, PSDCone())` or `@constraint(model, X <= 0, HermmitianPSDCone())`.
Stacktrace:
[...]

julia> @constraint(model, x >= y, Nonnegatives())
ERROR: At REPL[34]:1: `@constraint(model, x >= y, Nonnegatives())`: Unrecognized constraint building format. Tried to invoke `build_constraint(error, AffExpr[x[1,1] - 1 x[1,2] - 2; x[2,1] - 2 x[2,2] - 3], MathOptInterface.GreaterThan{Bool}(false), Nonnegatives())`, but no such method exists. This is due to specifying an unrecognized function, constraint set, and/or extra positional/keyword arguments.
[...]

PR3766

This PR will enable:

julia> model = Model();

julia> @variable(model, x[1:2, 1:2])
2×2 Matrix{VariableRef}:
 x[1,1]  x[1,2]
 x[2,1]  x[2,2]

julia> y = [1 2; 2 3]
2×2 Matrix{Int64}:
 1  2
 2  3

julia> @constraint(model, x >= y, Nonnegatives())
[x[1,1] - 1  x[1,2] - 2
 x[2,1] - 2  x[2,2] - 3] ∈ Nonnegatives()

julia> @constraint(model, x >= y, Nonpositives())
[x[1,1] - 1  x[1,2] - 2
 x[2,1] - 2  x[2,2] - 3] ∈ Nonpositives()

But it does not support:

julia> @constraint(model, x <= y, Nonnegatives())
ERROR: At REPL[7]:1: `@constraint(model, x <= y, Nonnegatives())`: The syntax `x <= y, Nonnegatives()` not supported. Use `y >= x, Nonnegatives()` instead.
Stacktrace:
[...]

This is because the user probably wants x - y in Nonpositives(), but it seems
ripe for confusion for the user to need to write Nonnegatives to get a
Nonpositives constraint.

Alternatively, like the current PSDCone support, we can make it
y - x in Nonegatives(), except then the primal and dual value will be the
negative of what the user "expects" (because the function is y - x not
x - y). There wasn't really a risk of this with PSDCone because people
really ever only mean the PSDCone. (We don't support a negative semidefinite
cone NSDCone.)

Your ask

Let's leave aside the question over the symbols ⪰ and ≥.

You're asking for

julia> @constraint(model, x >= y, ElementwiseComparison())

to be the same as @constraint(model, x >= y, Nonnegatives()) and

julia> @constraint(model, x <= y, ElementwiseComparison())

to mean x - y in Nonpositives() (i.e., the case we don't support).

Adding ElementwiseComparison is technically achievable, but it introduces
yet-another set into JuMP. At that point, you might as well just do x .<= y.

Tradeoffs

The current design is bad for someone who wants to write:

@constraint(model, x <= y)

to mean x - y in Nonpositives().

Their current alternatives are:

@constraint(model, x .<= y)
@constraint(model, y >= x, Nonnegatives())

We should probably also add:

@constraint(model, x - y in Nonpositives())  # TBD

Edit: oh, we can't add this because then the default x <= y would work...

Alternatives

Perhaps this issue is really just asking for an efficient
@constraint(model, x .<= y) for symmetric matrices.

Then there is no confusion anywhere.

If we do this, we should close JuliaLang/julia#3766 without merging.

I'll take a look at the details.

Confusions

I'm getting increasingly more and more confused.
so i'm not following what PSDCone has to do with it.

You suggested the example of erroring if >= is used with PSDCone. We can't do
that.

something else less friendly

Please assume that everyone in the JuMP community engages with the best positive
intentions. Asynchronous communication via short(-ish) comments is not always
a good medium for having subtle technical discussions... I think @matbesancon
holds the record for the most comments without a resolution:
jump-dev/MathOptInterface.jl#1023

[odow]

Okay. I remember why x .<= y was hard for Symmetric: it doesn't preserve under broadcasting

julia> using JuMP, LinearAlgebra

julia> model = Model();

julia> @variable(model, x[1:2, 1:2], Symmetric)
2×2 Symmetric{VariableRef, Matrix{VariableRef}}:
 x[1,1]  x[1,2]
 x[1,2]  x[2,2]

julia> x .- 1
2×2 Matrix{AffExpr}:
 x[1,1] - 1  x[1,2] - 1
 x[1,2] - 1  x[2,2] - 1

julia> y = LinearAlgebra.Symmetric([1 2; 2 3])
2×2 Symmetric{Int64, Matrix{Int64}}:
 1  2
 2  3

julia> y .- 1
2×2 Matrix{Int64}:
 0  1
 1  2

[LebedevRI]

Let's make this discussion more concrete.

Thank you.

What we currently support

<...>

This does match my understanding of the preexisting situation.

PR3766

This PR will enable:

julia> model = Model();

julia> @variable(model, x[1:2, 1:2])
2×2 Matrix{VariableRef}:
 x[1,1]  x[1,2]
 x[2,1]  x[2,2]

julia> y = [1 2; 2 3]
2×2 Matrix{Int64}:
 1  2
 2  3

julia> @constraint(model, x >= y, Nonnegatives())
[x[1,1] - 1  x[1,2] - 2
 x[2,1] - 2  x[2,2] - 3] ∈ Nonnegatives()

julia> @constraint(model, x >= y, Nonpositives())
[x[1,1] - 1  x[1,2] - 2
 x[2,1] - 2  x[2,2] - 3] ∈ Nonpositives()

Makes sense to me. It's an improvement on it's own already i think.

But it does not support:

julia> @constraint(model, x <= y, Nonnegatives())
ERROR: At REPL[7]:1: `@constraint(model, x <= y, Nonnegatives())`: The syntax `x <= y, Nonnegatives()` not supported. Use `y >= x, Nonnegatives()` instead.
Stacktrace:
[...]

This is because the user probably wants x - y in Nonpositives(), but it seems ripe for confusion for the user to need to write Nonnegatives to get a Nonpositives constraint.

YES! Extra cognitive load + confusion + much larger room for error, compared to the the non-matrix cases.

Alternatively, like the current PSDCone support, we can make it y - x in Nonegatives(), except then the primal and dual value will be the negative of what the user "expects" (because the function is y - x not x - y). There wasn't really a risk of this with PSDCone because people really ever only mean the PSDCone. (We don't support a negative semidefinite cone NSDCone.)

No comment.

Your ask

Let's leave aside the question over the symbols ⪰ and ≥.

You're asking for

julia> @constraint(model, x >= y, ElementwiseComparison())

to be the same as @constraint(model, x >= y, Nonnegatives()) and

julia> @constraint(model, x <= y, ElementwiseComparison())

to mean x - y in Nonpositives() (i.e., the case we don't support).

In other words, >= in @constraint(model, x >= y, Set()) is really more of a short-cut
to write @constraint(model, x-y in Set()), and <= symbol is not defined to do anything for matrices.
So what i'm essentially asking, is, when the set is ElementwiseComparison,
to treat >= as a comparison operator, and support it's inverse friend <=.
So i agree that your description matches my request.

Adding ElementwiseComparison is technically achievable, but it introduces yet-another set into JuMP.

Aha!

But it is still a bit not clear to me what is the problem with that,
what is the point of contention? Just the fact that said extra set is needed?

I haven't tried, but i would expect that such a set would require a very few changes,
and would basically immediately lower into the appropriate signedness set.

At that point, you might as well just do x .<= y.

Tradeoffs

The current design is bad for someone who wants to write:

@constraint(model, x <= y)

to mean x - y in Nonpositives().

Their current alternatives are:

@constraint(model, x .<= y)
@constraint(model, y >= x, Nonnegatives())

We should probably also add:

@constraint(model, x - y in Nonpositives())  # TBD

Yes.

Alternatives

Perhaps this issue is really just asking for an efficient @constraint(model, x .<= y) for symmetric matrices.

It does not, not really. Though such an improvement could be interesting.

There is an edge case though: what about vec(Symmetric(x_mat)) ?= vec(Symmetric(y_mat)),
can the duplicate constraints be avoided there? I'm guessing not?

Then there is no confusion anywhere.

If we do this, we should close JuliaLang/julia#3766 without merging.

I'll take a look at the details.

I was, in reality, thinking about (not necessarily symmetric) matrix support
for comparison operators, not for an efficient @constraint(model, x .<= y) for symmetric matrices.

Confusions

I'm getting increasingly more and more confused.
so i'm not following what PSDCone has to do with it.

You suggested the example of erroring if >= is used with PSDCone. We can't do that.

Err, that's not exactly correct. I did suggest that, but iff \succeq is introduced,
and iff it becomes the recommended way to explicitly specify a set,
so that >= no longer means \succeq, and >= becomes available
to be used as a comparison operator.
Since that is not happening, PSDCone is not as relevant.

something else less friendly

Please assume that everyone in the JuMP community engages with the best positive intentions. Asynchronous communication via short(-ish) comments is not always a good medium for having subtle technical discussions... I think @matbesancon holds the record for the most comments without a resolution: jump-dev/MathOptInterface.jl#1023

:) Right. I'm just saying that it was really NOT obvious to me if we were even understanding each other,
or talking past each other. The conversation became rather convoluted and very loopy/repetitive.

[odow]

I was, in reality, thinking about (not necessarily symmetric) matrix support for comparison operators

How is the existing x .<= y not sufficient? (Except for the redundancy with Symmetric and Hermitian)

[LebedevRI]

I suppose maybe it's not, other than stylistic preferences and differences in
printouts (which may be better for the .>= syntax anyways,
since it seems you can't print part of vector constraint).
But i may be forgetting something.

As i've said, it would be nice if matrix symmetry would be better respected in different constraint syntaxes,
but that seems rather not easier to implement than the subject in question (matrix inequalities).
One caveat would be, how will it be conveyed that symmetry has been exploited?

[odow]

I have added this to the agenda of this week's monthly developer call.

[LebedevRI]

@odow thank you! I guess i don't really care how the JuliaLang/julia#3765 is resolved, just that it would be nice if it is resolved.
I won't repeat the arguments again, maybe, indeed, making .[?]= work would be good-enough.

[blegat]

Yes, I agree making .== work would be a good solution

[odow]

Developer call says:

• We are not going to add support for ⪰ because the consensus was that majority of JuMP users will not perceive the difference, even if, yes, there is a technical difference that some users might expect.
• We are going to merge triangular dispatch (left-to-right template parameter chaining) JuliaLang/julia#3766
• We are not going to support x <= y, set because of the "negative dual" issue. The current support for PSDCone is regarded as a bug, but we're not going to fix because of the risk of breaking existing code.
• The .>= syntax is tricky, if not impossible to get working because of Julia's underlying broadcast behavior with symmetric matrices. We might add this in future if Base changes, but it isn't a priority.
• The recommendation for adding inequalities to symmetric matrices will be x >= l, Nonnegatives() and u >= x, Nonnegatives() which will enforce the correct (symmetric) element wise inequalities.

Closing this issue as won't-fix.

I'll leave JuliaLang/julia#3765 open until we merge JuliaLang/julia#3766.

[LebedevRI]

Thanks.

• We are not going to support x <= y, set because of the "negative dual" issue. The current support for PSDCone is regarded as a bug, but we're not going to fix because of the risk of breaking existing code.

Does this take into the account all of the lengthy ElementwiseComparison() set disscussion?
Given that .>= is effectively not happening either, that outcome is a bit surprising,
and all this UX disscussion was pointless.

[blegat]

We already have element-wise comparison with the .== syntax and we have conic comparison so #3766 seems natural and intuitive. Adding ElementwiseComparison() to work around the fact that broadcasting does not work for Symmetric matrices seems convoluted. What will we do once broadcasting is fixed for symmetric matrices in LinearAlgebra ? We don't want to add a feature that becomes useless or obsolete.

[LebedevRI]

What will we do once broadcasting is fixed for symmetric matrices in LinearAlgebra ?

Is there an upstream issue tracking that bug? #3770 (comment) made it sound like it's improbable that will be fixed.

[LebedevRI]

What will we do once broadcasting is fixed for symmetric matrices in LinearAlgebra ?

Is there an upstream issue tracking that bug?

JuliaLang/LinearAlgebra.jl#773 / JuliaLang/julia#19228 seem somewhat relevant,
but i'm not sure if that is the issue you've had in mind? @blegat

#3770 (comment) made it sound like it's improbable that will be fixed.

[LebedevRI]

I was, in reality, thinking about (not necessarily symmetric) matrix support for comparison operators

How is the existing x .<= y not sufficient? (Except for the redundancy with Symmetric and Hermitian)

Just for completeness, copying from jump-dev/MutableArithmetics.jl#295 (comment):
Apparently, because broadcasting doesn't just scalarize an "already-valid" vectorized operation,
but it may also broadcast whole dimensions to make the sizes match.
While this is nice, i hadn't even fully realized broadcast was doing that...

[odow]

Given A = (2, 3) and B = (2, 1, 3), the rules for broadcast are:

• Line up the two arrays

  A = (2, 3,)
  B = (2, 1, 3)
  

• pad any right dimensions to get same size

  A = (2, 3, 1)
  B = (2, 1, 3)
  

• repeat any dimensions to obtain a common size, so the last part of A gets repeated three times to make 2, 3, 3 and the second dimension of B gets repeated three times to make 2, 3, 3
• do the operation elementwise on arrays which are now the same size

This is notably NOT what you probably expected to happen, which is for broadcast to ignore the singleton dimension and return a size (2, 3) array.