Opened 9 months ago

Last modified 3 weeks ago

#29169 needs_review enhancement

Maximal and critical angles for cones

Reported by: mjo Owned by:
Priority: major Milestone: sage-9.2
Component: geometry Keywords:
Cc: jipilab, novoselt Merged in:
Authors: Michael Orlitzky Reviewers:
Report Upstream: N/A Work issues:
Branch: u/mjo/ticket/29169 (Commits) Commit: be2d48ec8a2a112e80178e00cb8614531365029b
Dependencies: #30605, #30654 Stopgaps:

Description

My paper with two algorithms for finding the maximal/critical angles between polyhedral convex cones was accepted:

http://www.optimization-online.org/DB_HTML/2019/01/7048.html

I've already written these algorithms in sage, so it would be nice to add them to the library. I just need to clean up some of the documentation.

Change History (39)

comment:1 Changed 9 months ago by mjo

Work in progress on u/mjo/ticket/29169. Doesn't build or anything yet.

comment:2 Changed 8 months ago by jipilab

  • Cc jipilab added

comment:3 Changed 8 months ago by mjo

  • Dependencies set to 26623

comment:4 Changed 6 months ago by mkoeppe

  • Milestone changed from sage-9.1 to sage-9.2

Batch modifying tickets that will likely not be ready for 9.1, based on a review of the ticket title, branch/review status, and last modification date.

comment:5 Changed 6 months ago by jipilab

  • Branch set to u/mjo/ticket/29169
  • Commit set to 99d77f01268d0d9e381956ba8f0a228a6ab4c19f

New commits:

c9437f1Trac #26623: add global entries for pending cone catalog citations.
792ac2bTrac #26623: add cones.<whatever> shortcuts for common cones.
ddf922cTrac #26623: add sage.geometry.cone_catalog to the documentation index.
f4a4d34Trac #26623: update sage.geometry.cone doctests to use the cone catalog.
9c3bf30Trac #29169: add the Or2020 reference.
6b73284Trac #29169: add the geometry/cone_critical_angles module.
08dd8aaTrac #29169: add the critical_angles and max_angle methods.
99d77f0UNDO ME: add the critical angles module to the doc index.

comment:6 Changed 4 months ago by git

  • Commit changed from 99d77f01268d0d9e381956ba8f0a228a6ab4c19f to fc255ba01aa5ce3d2f176aa3f217600c76b7a199

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

1edea71Trac #29169: add the Or2020 reference.
a2a025fTrac #29169: add the geometry/cone_critical_angles module.
4d15305Trac #29169: add the critical_angles and max_angle methods.
fc255baUNDO ME: add the critical angles module to the doc index.

comment:7 Changed 4 months ago by git

  • Commit changed from fc255ba01aa5ce3d2f176aa3f217600c76b7a199 to 016564a88267480c9d5f0a684845053a6cfc2576

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

468103bTrac #29169: add the Or2020 reference.
0481039Trac #29169: add the geometry/cone_critical_angles module.
b5ac943Trac #29169: add the critical_angles and max_angle methods.
a3c4535Trac #29169: fix a broken random_cone() test.
016564aUNDO ME: add the critical angles module to the doc index.

comment:8 Changed 4 months ago by mjo

  • Status changed from new to needs_review

I can keep tweaking the documentation forever, but I think this is in good condition to get some third-party feedback.

comment:9 Changed 4 months ago by mjo

(You can check the docs with sage -docbuild --include-tests-blocks --underscore)

comment:10 Changed 4 months ago by chapoton

  • Dependencies changed from 26623 to #26623
  • dependency field was missing #
  • Author field is empty

comment:11 Changed 4 months ago by mjo

  • Authors set to Michael Orlitzky
  • Dependencies #26623 deleted

comment:12 Changed 4 months ago by chapoton

patchbot complains about coverage and unused import

comment:13 Changed 4 months ago by git

  • Commit changed from 016564a88267480c9d5f0a684845053a6cfc2576 to 9e917d3730c1783f9f4507039420383284e394e6

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

468103bTrac #29169: add the Or2020 reference.
7d57406Trac #29169: add the geometry/cone_critical_angles module.
9c4af89Trac #29169: add the critical_angles and max_angle methods.
e4683f4Trac #29169: fix a broken random_cone() test.
9e917d3UNDO ME: add the critical angles module to the doc index.

comment:14 Changed 4 months ago by mjo

The unused import was legit (and is now fixed), but the coverage warning is for the critical_angles and max_angle functions that are wrapped by the methods of the same name in the cone class. All of the documentation (and unit tests) is on the cone methods, because that's what people will use.

(I can add some trivial tests if a happy coverage report is worth it, but so far I have avoided it as it's a complete waste of CPU time.)

comment:15 Changed 4 months ago by chapoton

I think that coverage is mandatory, sorry. Very simple tests would do, and you can tag them with #indirect doctest if they do not use the method directly.

More annoying, the latest patchbot reports a Timed Out, which is not acceptable.

comment:16 Changed 4 months ago by git

  • Commit changed from 9e917d3730c1783f9f4507039420383284e394e6 to 017cf1ad3a571dfe69c566b3b73057f2931644f6

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

4037b4eTrac #29169: add the geometry/cone_critical_angles module.
dfd587eTrac #29169: add the critical_angles and max_angle methods.
90a1914Trac #29169: fix a broken random_cone() test.
017cf1aUNDO ME: add the critical angles module to the doc index.

comment:17 Changed 4 months ago by mjo

I think that will make patchbot happy. I added some trivial tests where they were absent, and fixed the timeout issue (fingers crossed) by adding # long time in a bunch of places and reducing the max test size in others.

comment:18 Changed 6 weeks ago by gh-kliem

I'm willing to review this, but I'd much rather go through it one by one. Is it possible to implement the functions one by one?

I didn't check, maybe it all depends on each other.

comment:19 Changed 6 weeks ago by mjo

They're largely intertwined, but I could still commit some of the helper functions individually if that would make reviewing easier.

The reason for all of the separate functions (and the size of the module) is that critical_angles and max_angle do a lot of the same stuff. The max angle is a critical angle -- so critical_angles will find it -- but if all we're looking for is the biggest one then we can bail out as soon as we find it. That makes max_angle faster, but naturally the two share a lot of code, all of which I've tried to factor out into sensibly-named helper functions.

Once I defined the helper functions, they needed docstrings, tests, and examples.... meaning that a lot of relatively trivial things shared between critical_angles and max_angle take up a page or two of code. The inverse_spd function, for example, is two lines of implementation and a page of docs. If only one of critical_angles or max_angle were present, all of those functions would be inlined. Some of them make sense out of context, but others like check_gevp_feasibility are 100% tied to the implementation of critical_angles and max_angle, and are separate only to avoid copy/pasting their code.

comment:20 Changed 6 weeks ago by gh-kliem

Ok. I think I can also just look at the helper functions and work my way up.

comment:21 Changed 5 weeks ago by gh-kliem

  • Totally minor: If your in this situation, your cone is probably the trivial cone or the ambient space for real. It sounds like you claim that it isn't:
-            raise ValueError("this cone cannot be trivial or the ambient space")
+            raise ValueError("the cone may not be trivial or the ambient space")

or just delete the "this"

Other than this, the changes to src/sage/geometry/cone.py seem fine to me.

comment:22 Changed 5 weeks ago by mjo

I see what you're saying, but in informal English I think "may not" suffers from the same problem. Likewise for "must not". How about "the cone should not be..."? I have been unable to trick my brain into reading that incorrectly =)

comment:23 Changed 5 weeks ago by gh-kliem

the cone should not be is fine. I wouldn't want to worry about it too much.

comment:24 follow-up: Changed 5 weeks ago by gh-kliem

  • inverse_spd: Wouldn't that fit better into matrix/matrix2.pyx. This is where cholesky is defined.
  • _normalize_gevp_solution: Fine.
  • This error message is a full sentence.
    +        ValueError: The ambient dimension must be a positive integer (there
    +        are no nontrivial cones in dimension 0).
    
    I think just the part in the parenthesis would be a fine error message.

otherwise _random_admissible_cone is fine

  • +    # The standard cone containment tests don't work with normalized
    +    # vectors in QQbar.
    +    u_is_in_P = all( g.inner_product(u) >= 0 for g in P.dual() )
    +    v_is_in_Q = all( h.inner_product(v) >= 0 for h in Q.dual() )
    
    I would propose fixing this.
    sage: 3/2*vector([1,-1,0]) in P                                                                                                                                                     
    True
    sage: QQbar(2).sqrt()*vector([1,-1,0]) in P  
    ... # very long traceback
    ValueError: Cannot coerce non-integral Algebraic Real 1.414213562373095? to Integer
    
    This is awful and inconsistent. There is no way you can convert the point [3/2, -3/2, 0] to be a lattice point. According to the documentation it should have returned False. I would propose fixing this as a dependency of this ticket:
    • outsource dual_rays or similar from dual
    • use dual_rays to improve _contains to behave coercion-like with your code
    • then you can just compute u in P.

If I ask whether a rational or real algebraic point is contained in a cone, I want coercion.

The second thing: all( z.inner_product(h) >= 0 for h in Q ) could then be reduced to z in Q.dual() (which is a bit wasteful, but not more than the current state). This could be exposed somehow as well as this is the general question whether Q satisfies some inequality, but that is probably not a priority.

  • gevp_licis not fast, but if the time of running this is not relevant, it should be fine.

comment:25 Changed 5 weeks ago by gh-kliem

In examples you define this G_index_sets = list(gevp_licis(G)) without ever using it or printing it.

comment:26 in reply to: ↑ 24 ; follow-up: Changed 5 weeks ago by mjo

  • Cc novoselt added

Replying to gh-kliem:

This is awful and inconsistent... If I ask whether a rational or real algebraic point is contained in a cone, I want coercion.

Me too, but I anticipate that Andrey will object. The cone module is fairly strict; if you want to do anything "dirty," you have to leave the lattice and move to your own vector space. You've probably seen this in other tickets of mine where I have to throw away the lattice structure of (say) a dual cone to be able to treat its rays as objects in the primal space. It's often annoying, but it's consistent within the cone module. And it's technically correct, even from the convex analysis/geometry perspective, so long as you use the traditional definition of "dual cone" where the elements lie in a dual vector space. So I've learned to live with it.

In either case, you are right that it should return False as documented.

comment:27 Changed 5 weeks ago by git

  • Commit changed from 017cf1ad3a571dfe69c566b3b73057f2931644f6 to 21f1cd2b9ce9838adb7f3f8cda6258dd62e407c2

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

c9cb3d5Trac #29169: add the Or2020 reference.
583967dTrac #29169: add the geometry/cone_critical_angles module.
61a46dcTrac #29169: add the critical_angles and max_angle methods.
c323228Trac #29169: fix a broken random_cone() test.
ec72863UNDO ME: add the critical angles module to the doc index.
f32eb63Trac #29169: tweak two max/critical angle error messages.
c8dd547Trac #29169: improve error message for _random_admissible_cone().
21f1cd2Trac #29169: remove two unused example lines.

comment:28 in reply to: ↑ 26 ; follow-up: Changed 5 weeks ago by gh-kliem

How about def contains(self, coerce=False) or similar. So Andrey can have his way, but we are still not creating copies of stuff at places where they do not belong.

It's still strange that it would accept a vector but not a point in the dual lattice as input.

In a way it makes sense:

sage: 2 in range(2,5)                                                                                                                                                                                                                                                                                                                                                      
True
sage: 2.5 in range(2,5)                                                                                                                                                                                                                                                                                                                                                    
False

This is completely natural and I can see that you want this behavior for a cone. However, the other question should still be accessible. As you cannot change the base ring of a cone (can you?) it would be nice for the cone to provide some answer for you.

Replying to mjo:

Replying to gh-kliem:

This is awful and inconsistent... If I ask whether a rational or real algebraic point is contained in a cone, I want coercion.

Me too, but I anticipate that Andrey will object. The cone module is fairly strict; if you want to do anything "dirty," you have to leave the lattice and move to your own vector space. You've probably seen this in other tickets of mine where I have to throw away the lattice structure of (say) a dual cone to be able to treat its rays as objects in the primal space. It's often annoying, but it's consistent within the cone module. And it's technically correct, even from the convex analysis/geometry perspective, so long as you use the traditional definition of "dual cone" where the elements lie in a dual vector space. So I've learned to live with it.

In either case, you are right that it should return False as documented.

comment:29 in reply to: ↑ 28 ; follow-up: Changed 5 weeks ago by novoselt

Replying to gh-kliem:

It's still strange that it would accept a vector but not a point in the dual lattice as input.

For any vector space with a basis it is quite natural to have conversion from the set of numbers, i.e. a vector from ZZ^n to that vector space. And the backward map is fairly natural as well. However a composition of such maps does not give a natural map between two vector spaces with chosen bases. On that ground dual lattice elements are not accepted as points in the primal lattice. And the reason why this particular case is explicitly prohibited is that it is a fairly common mistake in toric geometry to work with elements/cones of the wrong lattice, so it helps to catch mistakes. Volker's alternative code was actually having completely unrelated classes for a lattice and its dual to make sure that there is no conversion/coercion happening between them.

So I am strongly against changing this behaviour in any situation. Having explicit non-default way to mix lattices is perfectly fine, of course.

comment:30 follow-up: Changed 5 weeks ago by novoselt

Regarding checking rational/algebraic/float vectors for cone containment - it would be nice for it to work. The problem is, I guess, how to write down the check. For my purposes working with rational numbers was always sufficient and it has the advantage of exact computations, so no need to worry about precision issues. If you can add code to handle more cases - please do it and it has nothing to do with distinction between dual lattices!

comment:31 in reply to: ↑ 29 ; follow-up: Changed 5 weeks ago by mkoeppe

Replying to novoselt:

[...] dual lattice elements are not accepted as points in the primal lattice. And the reason why this particular case is explicitly prohibited is that it is a fairly common mistake in toric geometry to work with elements/cones of the wrong lattice, so it helps to catch mistakes. Volker's alternative code was actually having completely unrelated classes for a lattice and its dual to make sure that there is no conversion/coercion happening between them.

As a side remark, a weakness of the vector space / free module code in sage.modules is that dual spaces / modules are not implemented. The newer code in sage.tensor.modules fixes this, and it would be a great improvement if sage.geometry.cone and sage.geometry.toric_lattice could integrate with that.

comment:32 in reply to: ↑ 30 Changed 5 weeks ago by mjo

Replying to novoselt:

If you can add code to handle more cases - please do it and it has nothing to do with distinction between dual lattices!

I was mainly worried about being able to take, say, M(1,1) and (through coercion) come to the conclusion that it lives in the nonnegative quadrant of the lattice N. But thanks to some implementation details, that appears not to be a problem. I'll post something soon.

comment:33 Changed 5 weeks ago by mjo

  • Dependencies set to #30605

comment:34 in reply to: ↑ 31 Changed 5 weeks ago by novoselt

Replying to mkoeppe:

As a side remark, a weakness of the vector space / free module code in sage.modules is that dual spaces / modules are not implemented. The newer code in sage.tensor.modules fixes this, and it would be a great improvement if sage.geometry.cone and sage.geometry.toric_lattice could integrate with that.

Thank you so much for the pointer! I surely wished for this when I was implementing toric stuff.

comment:35 Changed 5 weeks ago by git

  • Commit changed from 21f1cd2b9ce9838adb7f3f8cda6258dd62e407c2 to adeb746e1e5628bc9bb8d9d548dd70889a881069

Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:

f8b8b37Trac #30605: update polyhedral cone containment documentation.
db132e0Trac #29169: add the Or2020 reference.
3ff8035Trac #29169: add the geometry/cone_critical_angles module.
68f339bTrac #29169: add the critical_angles and max_angle methods.
358aec4Trac #29169: fix a broken random_cone() test.
fd5ac83UNDO ME: add the critical angles module to the doc index.
f4fc161Trac #29169: tweak two max/critical angle error messages.
64c0784Trac #29169: improve error message for _random_admissible_cone().
31b47afTrac #29169: remove two unused example lines.
adeb746Trac #29169: take advantage of improvements in trac #30605.

comment:36 Changed 4 weeks ago by mjo

  • Dependencies changed from #30605 to #30605, #30654

comment:37 Changed 4 weeks ago by mjo

  • Status changed from needs_review to needs_work

I'll need to rebase this again once #30654 is reviewed.

comment:38 Changed 3 weeks ago by git

  • Commit changed from adeb746e1e5628bc9bb8d9d548dd70889a881069 to be2d48ec8a2a112e80178e00cb8614531365029b

Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:

aa506fcTrac #30654: mention inverse_positive_definite() in inverse().
5a834d2Trac #29169: add the Or2020 reference.
4341a39Trac #29169: add the geometry/cone_critical_angles module.
cfce9c1Trac #29169: add the critical_angles and max_angle methods.
970db5dTrac #29169: fix a broken random_cone() test.
714e3b7UNDO ME: add the critical angles module to the doc index.
f72376eTrac #29169: tweak two max/critical angle error messages.
511c360Trac #29169: improve error message for _random_admissible_cone().
4b452d5Trac #29169: remove two unused example lines.
be2d48eTrac #29169: take advantage of improvements in trac #30605.

comment:39 Changed 3 weeks ago by mjo

  • Status changed from needs_work to needs_review

Cherry-picked the two commits from #30654, and replaced the inverse_spd function with the inverse_positive_definite method.

Note: See TracTickets for help on using tickets.