[WIP] Evaluation distance to set #1023

matbesancon · 2020-02-11T10:44:36Z

Fairly straightforward, just like eval_variables evaluates the value of a function with a given value mapping of the variables, this PR adds the evaluation of values belonging to a set.

I reused the Base.in function, which seemed appropriate here since it's not used elsewhere with this signature.

One thing to work out as pointed by @mtanneau is the propagation of absolute and relative tolerance

codecov-io · 2020-02-11T11:41:48Z

Codecov Report

Merging #1023 into master will decrease coverage by 0.19%.
The diff coverage is 91.3%.

@@           Coverage Diff            @@
##           master   #1023     +/-   ##
========================================
- Coverage   95.19%     95%   -0.2%     
========================================
  Files         100     100             
  Lines       11279   11544    +265     
========================================
+ Hits        10737   10967    +230     
- Misses        542     577     +35

Impacted Files	Coverage Δ
src/sets.jl	`80.84% <91.3%> (-14.45%)`	⬇️
src/Bridges/Constraint/geomean.jl	`95.23% <0%> (-3.78%)`	⬇️
src/FileFormats/MOF/nonlinear.jl	`100% <0%> (ø)`	⬆️
src/Bridges/Constraint/det.jl	`100% <0%> (ø)`	⬆️
src/Test/intlinear.jl	`100% <0%> (ø)`	⬆️
src/Test/contconic.jl	`99.29% <0%> (+0.02%)`	⬆️
src/Utilities/model.jl	`94.91% <0%> (+0.04%)`	⬆️
src/Utilities/mutable_arithmetics.jl	`90.09% <0%> (+0.09%)`	⬆️
src/FileFormats/MOF/write.jl	`66.31% <0%> (+0.72%)`	⬆️
... and 2 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 93cdc8a...b19f884. Read the comment docs.

mlubin · 2020-02-11T14:28:31Z

IMO for this to be generally useful it needs to support tolerances in some form.

matbesancon · 2020-02-11T14:39:00Z

Agreed, I wasn't sure yet how to propagate them, rtol and atol keywords similar to isapprox?

odow · 2020-02-11T14:39:43Z

I would prefer a distance function which returned the distance from a point to feasibility inside the set.

function distance(x::T, set::LessThan{T}) where {T}
    return max(zero(T), x - set.upper)
end

function distance(x::T, set::Integer) where {T}
    return abs(x - round(Int, x))
end

function distance(x::Vector{T}, set::SOS1{T}) where {T}
    _, index = findmax(abs.(x))
    return sqrt(sum(x[i]^2 for i = 1:length(x) if i != index))
end

Feasibility checkers could then use these distances and some notion of tolerances.

matbesancon · 2020-02-11T14:42:10Z

a notion of distance makes sense for convex sets, integer sets get more ambiguous

odow · 2020-02-11T14:55:48Z

integer sets get more ambiguous

Sure. But if we have a notion of belonging to a set with a tolerance, then we need some corresponding distance to check against.

For some sets, e.g., SOSII, it isn't immediately obvious what a good distance is. L2 norm of the elements excluding the two largest adjacent absolute values?

We might even consider returning a vector to feasibility, rather than the norm. So a PSD constraint is the vector of eigenvalues min.(0.0, eigvals(X)).

That would let the feasibility checker choose an appropriate norm as well. In the SOSII case, it isn't apparent if you want to use the L2 norm.

mtanneau · 2020-02-11T15:00:11Z

integer sets get more ambiguous

That can be the distance from x to the lattice, rather than from x to a continuous set.

Then, as odow seems to suggest, which norm to use could be up to the user.
In a flavor similar to norm(x, 2) or norm(x, Inf).

mtanneau · 2020-02-11T15:16:47Z

One (probably obvious) point, maybe useful for the docs: for a constraint f(x) ∈ S, here we are talking of "feasibility" as the distance from f(x) to S.
The distance from x to the feasible set X = {x | f(x) ∈ S} can be something else entirely.

Thus, here we are really addressing the question of "how close is this constraint to being satisfied". (which is 100% fine by me)

matbesancon · 2020-02-11T15:28:51Z

just to be sure: you mean we deal with distance in the image space, and not in the input space right?

matbesancon · 2020-02-11T15:33:17Z

One more corner case, in the case of dimension mismatch, I simply returned false for set belonging. When using a distance, this would be +Inf I guess

matbesancon · 2020-02-11T15:35:26Z

We might even consider returning a vector to feasibility, rather than the norm. So a PSD constraint is the vector of eigenvalues min.(0.0, eigvals(X)).

+1

mtanneau · 2020-02-11T15:37:01Z

in the case of dimension mismatch

Why not throw a DimensionMismatch?

matbesancon · 2020-02-11T15:49:50Z

Not a fan of proliferating thrown errors all over the place, but yeah seems appropriate with this change

matbesancon · 2020-02-11T16:06:23Z

I've started a re-write up to NormSpectralCone.
One open question is whether to use min(a, b, f(x) - S) when requiring, say, a >= 0, b >= 0, f(x) in S

matbesancon · 2020-02-11T16:09:18Z

note here that _distance returns an unsigned distance, which is then converted in set_distance

src/Utilities/set_eval.jl

odow

I would prefer that we define a distance_from_set(x, set) function and include them along-side the set definitions in set.jl. Then it's obvious which sets have distances and which don't.

If x is inside the set, the distance should be 0. We shouldn't measure distance from infeasibility.

src/Utilities/set_eval.jl

matbesancon · 2020-02-11T18:17:46Z

@odow the _distance function remains internal and computed the unsigned distance to avoid checking if >= 0 in every single method.
The external function is the top one, computing the distance, and then returning max(0, dist)

odow · 2020-02-11T18:53:21Z

I understand the current split. I'm suggesting a public method for each set, with no re-direction to the zero-checker. Maybe a few more lines of code, but easier to understand for each set in isolation.

matbesancon · 2020-02-11T19:36:20Z

Ok yes I see your point, will do that. I still think we should keep the non-checked version, which would be comparable to a noargcheck or @inbounds``

odow · 2020-02-11T19:40:23Z

The zero check is essentially free. No one should be using this in any performance-critical code, and even if they are, max(0, x) isn't going to be the bottleneck. I'm in favor of simple as possible, one-way-to-do-things. You should wait for someone else to chime in though, in case I'm by myself on this.

matbesancon · 2020-02-11T22:21:12Z

Maybe a few more lines of code, but easier to understand for each set in isolation.

the bother with this is mostly that it means maintaining boilerplate over all sets instead of just at one place. That's why one function could do the dimension check, possibly throwing an error, call the set-specific method and return the cleaned-up version, otherwise these three steps have to be copy-pasted all over the place

matbesancon · 2020-03-27T16:24:30Z

weird sets can define their distance in whatever "best" way for them (indicator constraints are an example)

matbesancon · 2020-05-16T15:05:45Z

Open questions for distances to set:
To extend the distance system, do we need a fallback for a new distance to the default one?
To extend sets, do you extend only the default distance?

odow · 2020-05-16T16:04:04Z

Yes, presumably we have a fallback like:

function distance_to_set(::AbstractDistanceFunction, x, s)
    return distance_to_set(DefaultDistance(), x, s)
end

Note that DistanceFunction should probably be renamed AbstractDistanceFunction.

…my test distance

matbesancon · 2020-05-16T16:07:48Z

Another issue comes from complementarity constraints:

MOI.Complements() requires access to variable bounds to determine if a given vector of values
respects the complementarity constraint, which is not possible with the current signature.
A possible fix would be accessing the variable bounds in a generic way.
Another would be to specify complementarity constraints with the associated bounds Complementarity(Vector{<:AbstractSet}).

odow · 2020-05-16T18:38:44Z

The main reason why complementary bounds are defined separately is to avoid issues like the following, where declare bounds twice.

model = Model()
@variable(model, 1 <= x <= 2)
@constraint(model, 2x - 1 ⟂ { x >= 0})

It seems reasonable that we wouldn't implement the DefaultDistance for complements. They're a bit if a special case. You could imagine defining a new distance function that required the l and u bound vectors as input.

matbesancon · 2020-05-16T20:17:17Z

Yes agreed, even though it makes Complements a bit of a corner case for many things outside of distance, because of this non-self-contained constraint

odow

I think we still need to nail down the theoretical concept of the distance we are measuring, particularly for some of the conic sets with strict inequalities.

For some sets you can get closed-form expressions for the minimum L2 distance (e.g., #1023 (comment)), and I think we should use these where possible.

odow · 2020-05-19T13:50:00Z

src/sets.jl

@@ -1,6 +1,8 @@
 # Sets

 # Note: When adding a new set, also add it to Utilities.Model.
+import LinearAlgebra
+using LinearAlgebra: dot


Why import dot but not norm2? Just prefix everything with LinearAlgebra.

odow · 2020-05-19T13:50:36Z

src/sets.jl

+    distance_to_set(v, s)
+
+Compute the distance of a value to a set.
+For some vector-valued sets, can return a vector of distances.


This is still an issue.

odow · 2020-05-19T13:51:10Z

src/sets.jl

+
+Distance function used to evaluate the distance from a point to a set.
+New subtypes of `AbstractDistance` must implement fallbacks for sets they don't cover and implement
+`distance_to_set(::Distance, v, s::S)` for sets they override the distance for.


I thought we were providing the AbstractDistance fallback?

This is still an issue.

mistake on my part this should be removed

I thought we were providing the AbstractDistance fallback?

we can provide one d(::AbstractDistance, v, s) = d(::DefaultDistance, v, s) but not be more specific because of method ambiguity

odow · 2020-05-19T13:56:16Z

src/sets.jl

+    t = v[1]
+    xs = v[2:end]
+    result = LinearAlgebra.norm2(xs) - t
+    return max(result, zero(result)) # avoids sqrt


What is the comment # avoids sqrt for?

I have a general unresolved concern over our choice of distance. It seems like this is the "epigraph violation" distance, which is not the same as the Euclidean distance to the set. What is the precedent for how other solvers measure constraint feasibility? (@chriscoey / @lkapelevich how does Hypatia?)

I personally use the epigraph violation d = x1 - norm(x[2:end], 2), mainly because

it's easy to compute

it tells by how much the conic constraint should be "relaxed" for the point to be feasible (where "relaxed" means that you shift the cone's apex)

if the SOC constraint is written as x2^2 + ... xn^2 <= x1^2, the epigraph violation is the square root of this quadratic constraint's absolute violation.

Looks like COSMO uses euclidean distance for computing projections.

What about sets like the exponential cone needing domain constraints like y > 0?

What is the comment # avoids sqrt for?

legacy comment, just removed. For domain constraints it can be an open question, because some parts of the distance may be non-computable / infinite or complex if the domain is not respected

Okay, so for conic sets "Epigraph violation or Inf if point fails feasibility bounds" seems like an easily computable and sensible definition.

And you call it EpigraphDistance or something like...

For norm two distances you can use what proximal solvers do:
https://web.stanford.edu/~boyd/papers/pdf/prox_algs.pdf page 183 (section 6.3)

So, SOC is easy, PSD you just need eigen decomposition, Power and Exponential are well posed even with the extra constraints, but need a Newton...

odow · 2020-05-19T13:57:22Z

src/sets.jl

+    u = v[2]
+    xs = v[3:end]
+    return LinearAlgebra.norm2(
+        (max(-t, zero(t)), max(-u, zero(u)), max(dot(xs,xs) - 2 * t * u))


I think we need to think this through a little. What space is this a distance in?

The core issue here is that the cone is characterized by three constraints.
With vector-distance, this would have been:

max(-u, 0) max(-v, 0) max(-2tu + norm_squared(x) , 0)

since we don't have that, the l2-norm of these three epigraph violations is the thing computed here.

odow · 2020-05-19T14:04:34Z

src/sets.jl

+            return d
+        end
+    end
+end


I don't think this function is correct. Something like:

function distance_to_set(::DefaultDistance, v::AbstractVector{T}, ::SOS1{T}) where {T <: Real} _check_dimension(v, s) _, i = findmax(abs.(v)) return LinearAlgebra.norm2([v[j] for j = 1:length(v) if j != i]) end

chriscoey · 2020-05-19T17:17:10Z

My comment here isn't actually very helpful but: Hypatia doesn’t use a notion of distance to set anywhere. It’s tough because there are many notions of distance and the right one depends on context. Also not all conic sets are typically defined by epigraph/hypographs, e.g. the PSD cone or the set of nonnegative polynomials or the copositive cone.

chriscoey · 2020-05-19T17:20:29Z

Also for many sets the Euclidean distance can be hard to compute (eg for exponential cone in general, though some "cases" are just violation on variable nonnegativity)

joaquimg · 2020-05-19T23:21:18Z

I will put it here, so it does not gets lost in the inner conversation:

For norm two distances you can use what proximal solvers do:
https://web.stanford.edu/~boyd/papers/pdf/prox_algs.pdf page 183 (section 6.3)

So, SOC is easy, PSD you just need eigen decomposition, Power and Exponential are well posed even with the extra constraints, but need a Newton...

joaquimg · 2020-05-19T23:28:18Z

Many of these projections are implemented here: https://github.com/kul-forbes/ProximalOperators.jl/tree/master/src/functions

odow · 2020-05-20T15:04:02Z

Considering this PR is approaching the longest PR in the history of MOI (which @matbesancon also holds first and second place in #709 and #877), maybe this is a sign it should be in a package first. Then @matbesancon can experiment with the different set definitions in peace, not worry about pulling in heavy dependencies for the proximal operators or the SVD, and see what works.

Edit: I will note, beginning the PR with "Fairly straightforward," was an understatement 😆

matbesancon · 2020-05-20T15:08:51Z

maybe this is a sign it should be in a package first

do you mean the distances themselves or the advanced features?

Dear 2-month-ago self, what what you done?

odow · 2020-05-20T15:21:28Z

I mean the distances, the feasibility checker, everything. Call it matbesancon/InfeasibilityToolkit.jl.

This approach worked for Dualization.jl, and we brought the dual_set definitions in once they stabilized.

matbesancon · 2020-05-20T15:25:37Z

I should have a skeleton more or less there for Friday's meeting

matbesancon · 2020-05-21T18:42:58Z

So, with the heart broken, we can close this PR which is replaced by https://github.com/matbesancon/MathOptSetDistances.jl

matbesancon added 2 commits February 11, 2020 10:55

first sets

b606dc3

more sets

6670920

matbesancon changed the title ~~Evaluation of values belonging to sets~~ [WIP] Evaluation of values belonging to sets Feb 11, 2020

forgot subslice

93ac7f9

rewrite as distance function

55ed9d6

sign confusion

09663dd

mtanneau reviewed Feb 11, 2020

View reviewed changes

src/Utilities/set_eval.jl Outdated Show resolved Hide resolved

added integers

9d4def8

odow reviewed Feb 11, 2020

View reviewed changes

src/Utilities/set_eval.jl Outdated Show resolved Hide resolved

src/Utilities/set_eval.jl Outdated Show resolved Hide resolved

src/Utilities/set_eval.jl Outdated Show resolved Hide resolved

src/Utilities/set_eval.jl Outdated Show resolved Hide resolved

function for verification not being toplevel

c02e1ba

matbesancon added 2 commits March 31, 2020 15:42

Merge branch 'master' into eval-set

ed1d33b

Merge branch 'master' into eval-set

f80f399

default distance

54ff584

abstract and default distance for sets, added complementarity and dum…

68888b5

…my test distance

remove complementarity

21c759f

increase tolerance

9a2e898

odow reviewed May 19, 2020

View reviewed changes

update doc

13725c3

matbesancon closed this May 21, 2020

odow mentioned this pull request Aug 21, 2020

Add value for constraints with a custom var_value function jump-dev/JuMP.jl#2317

Merged

matbesancon deleted the eval-set branch February 9, 2021 16:44

odow mentioned this pull request Feb 15, 2021

First pass at a feasibility checker. jump-dev/JuMP.jl#2466

Merged

odow mentioned this pull request Oct 23, 2022

Implement Utilities.distance_to_set for more sets #2033

Closed

odow mentioned this pull request Nov 14, 2022

[Utilities] add distance_to_set #2048

Merged

odow mentioned this pull request Jun 11, 2024

Set inequality syntax suggestion jump-dev/JuMP.jl#3770

Closed

[WIP] Evaluation distance to set #1023

[WIP] Evaluation distance to set #1023

Conversation

matbesancon commented Feb 11, 2020

codecov-io commented Feb 11, 2020 • edited

Codecov Report

mlubin commented Feb 11, 2020

matbesancon commented Feb 11, 2020

odow commented Feb 11, 2020 • edited

matbesancon commented Feb 11, 2020

odow commented Feb 11, 2020

mtanneau commented Feb 11, 2020

mtanneau commented Feb 11, 2020

matbesancon commented Feb 11, 2020

matbesancon commented Feb 11, 2020

matbesancon commented Feb 11, 2020

mtanneau commented Feb 11, 2020

matbesancon commented Feb 11, 2020

matbesancon commented Feb 11, 2020

matbesancon commented Feb 11, 2020

odow left a comment

Choose a reason for hiding this comment

matbesancon commented Feb 11, 2020

odow commented Feb 11, 2020

matbesancon commented Feb 11, 2020

odow commented Feb 11, 2020

matbesancon commented Feb 11, 2020

matbesancon commented Mar 27, 2020

matbesancon commented May 16, 2020

odow commented May 16, 2020

matbesancon commented May 16, 2020 • edited

odow commented May 16, 2020

matbesancon commented May 16, 2020

odow left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

chriscoey commented May 19, 2020

chriscoey commented May 19, 2020

joaquimg commented May 19, 2020

joaquimg commented May 19, 2020

odow commented May 20, 2020 • edited

matbesancon commented May 20, 2020

odow commented May 20, 2020

matbesancon commented May 20, 2020

matbesancon commented May 21, 2020

codecov-io commented Feb 11, 2020 •

edited

odow commented Feb 11, 2020 •

edited

matbesancon commented May 16, 2020 •

edited

odow commented May 20, 2020 •

edited