Constructions for common polyhedral cones

[orlitzky]

When working with graphs, we have a nice tab-completed list of common graphs that we can construct; for example,

sage: graphs.BuckyBall()
Bucky Ball: Graph on 60 vertices

More and more, I find myself wishing we had the same thing for convex cones, particularly in test cases. I have constructions for all of the following and use them regularly:

• The nonnegative orthant.
• The full ambient space.
• The schur cone that induces the majorization ordering.
• The permutation-invariant "rearrangement" cone.
• The trivial cone, since trivial_cone(n) looks a lot nicer than Cone([], ToricLattice(n)).

et cetera; I'm sure I can think of more. Each of these makes sense in any dimension n.

I'm wondering if you think this would be a useful feature to make available as e.g.

sage: cones.nonnegative_orthant(5)
5-d cone in 5-d lattice N

They should all probably take a "lattice" argument too now that I think about it.

CC: @novoselt

Component: geometry

Author: Michael Orlitzky

Branch/Commit: b2aab91

Reviewer: Jonathan Kliem

Issue created by migration from https://trac.sagemath.org/ticket/26623

[novoselt]

comment:1

They certainly should take a lattice argument! From a toric perspective it is probably more natural to have these as methods of lattices, i.e. something like

ToricLattice(n).trivial_cone()

although of course this makes it impossible to use ZZ^n. And perhaps even more natural is to have not cones, but associated toric varieties which are already available:

sage: toric_varieties.A(5).fan().generating_cone(0).rays()
N(1, 0, 0, 0, 0),
N(0, 1, 0, 0, 0),
N(0, 0, 1, 0, 0),
N(0, 0, 0, 1, 0),
N(0, 0, 0, 0, 1)
in 5-d lattice N

[orlitzky]

comment:2

Replying to @novoselt:

They certainly should take a lattice argument! From a toric perspective it is probably more natural to have these as methods of lattices, i.e. something like

ToricLattice(n).trivial_cone()

although of course this makes it impossible to use ZZ^n.

I think from a design standpoint, that's where the method belongs, but there are two things that bug me about it.

First: in many cases, the user won't care about the lattice and will want "the default." This is basically what happens when you do ToricLattice(n).nonnegative_orthant(), but it feels a little weird to specify a particular lattice when you don't care about it. By analogy, the Cone() function should also be a method on a ToricLattice object. A lot of the time you just don't care though, so it's nice to be able to get a sensible default without having to construct a generic lattice of the correct size first. On the other hand, I like that this eliminates the n and lattice parameters from the cone-creating functions. And it would ensure that (for example) if you only create one lattice object, then you can intersect any two cones created by methods on it. I'll have to think harder about this but it may be the way to go.

Second: tab completion on graphs.<tab> is so nice. Having a list of twenty or so cones interspersed in the tab-completion list for lattices would be bad for both people who want to see lattice methods and people who want to see predefined cones.

(Might we have something like lattice.cones.<tab> where lattice is some ToricLattice?)

And perhaps even more natural is to have not cones, but associated toric varieties which are already available:

sage: toric_varieties.A(5).fan().generating_cone(0).rays()
...

My main objection to this is that I'm not smart enough to implement it or remember which cone is associated with which toric variety =)

Having the code above be the implementation for nonnegative_orthant(5) would be fine with me though.

[novoselt]

comment:3

"A lot of the time you just don't care though" - I always care, we are just using cones for different things ;-) But you are right - just because cones are associated to lattices does not mean that they have to be methods of lattices. So I think cones.etc is the right interface, but there has to be an option to specify the lattice and the default should be the same as for Cone(...).

[orlitzky]

comment:4

I have another question. There are a few vector spaces associated with each cone, that can be thought of as cones themselves:

• The span of the cone K.
• Its orthogonal complement, (span(K))-perp.
• Its lineality space (the largest subspace it contains).

The span of the cone and its lineality space are already available as,

sage: K = Cone([(1,1),(0,1)])
sage: K.span()
Sublattice <N(1, 0), N(0, 1)>
sage: span(K)
Free module of degree 2 and rank 2 over Integer Ring
Echelon basis matrix:
[1 0]
[0 1]
sage: K.linear_subspace()
Vector space of degree 2 and dimension 0 over Rational Field
Basis matrix:
[]

But none of those return cones directly. Having a cone is handy when you want to compute e.g. the intersection of one of these things with some other cone.

If we were to add similar methods that return cones, how would you prefer to do it? We could either...

1. Add methods to the cones themselves, for things like K.span_cone(), K.lineality_space(), K.orthogonal_complement()
2. Add them as part of the new cones.<whatever> interface, with something like cones.span_of(K).

The first one seems more correct to me, but risks polluting the namespace with a bunch of methods that all look the same but return different things.

The second one does make it obvious that you're going to get a cone back, but only if people know to look there for what should in principle be a method on the cone itself.

[orlitzky]

comment:5

Option 3:

For the span and the orthogonal complement, we could sensibly add a cone() method somewhere in the module hierarchy, like

sage: K.orthogonal_sublattice().cone()

which would be a wrapper around

sage: Cone( c*g for g in K.orthogonal_sublattice().gens() for c in [-1, 1] )

The problem with that is, when starting from a cone and passing through a vector space, the lattice information is lost. If I start with a cone in a lattice M, take its linear_subspace(), and then turn it into cone...

sage: K = Cone([(1,0,0),(-1,0,0),(0,1,0)], ToricLattice(3,'M'))
sage: Cone(K.linear_subspace().gens())
1-d cone in 3-d lattice N

it should live in the same lattice as the original cone. I don't see how to fix that immediately.

[novoselt]

comment:6

Just a reminder - orthogonality in the toric setting is handled via dual lattices:

sage: c = Cone([(1,2,3), (4,5,6)])
sage: c
2-d cone in 3-d lattice N
sage: c.orthogonal_sublattice()
Sublattice <M(1, -2, 1)>
sage: c.span()
Sublattice <N(1, 2, 3), N(0, 3, 6)>

I can't imagine use cases where I would like to treat all these spaces as cones. Can you please explain it in more details, i.e. what is wrong with treating them as, well, spaces? If you really do want to have cones associated to spaces, then you can easily attach a method to toric lattices and sublattices, or, alternatively or in addition, you can modify the Cone constructor to handle lattices and spaces. I just tried Cone(c.span()) for the above example and probably hit some infinite recursion that is in addition very slow - after a few minutes it is still doing something.

But I am still unclear why would one want to turn a space into a cone, what is there to gain from such a representation?

[orlitzky]

comment:7

Oof, it looks like I'm missing trac emails again.

It would just be nice to have the cone methods available on subspaces. For example, Corollary 1 in these notes,

http://michael.orlitzky.com/notes/on_z-operators_and_viability_theorems.pdf

is very easy to test as the intersection of a closed convex (dual) cone and a subspace. If span(x) is a cone, and if its orthogonal complement is a cone, then I can use the cone intersection method that we already have to test it.

Another example is the "maximal angle between two convex cones." This is a newish concept (from 2016), and I've implemented it as a method on cones. However, the "principal angle between subspaces" is a classical idea going back to the 1950s, and they turn out to be equivalent for cones that are subspaces. If subspaces can be treated as cones, then I can use the same method that I have on cones to compute principal angles between subspaces. And so on.

The reason Cone(c.span()) "hangs" is because it carefully enumerates every element of the vector space.

I think there's plenty of work to do just to get the trivial cone, nonnegative orthant, and a few others implemented and documented. I'll focus on that and worry about the other questions later.

[orlitzky]

Commit: eb072ce

[orlitzky]

comment:8

Here's my first attempt at this. I'm open to other suggestions if you don't like the module or class names.

---

New commits:

6645911	Trac #26623: add global entries for pending sage.geometry.cones citations.
84f5635	Trac #26623: add cones. shortcuts for common cones.
dd759dc	Trac #26623: add the sage.geometry.cones module to the documentation index.
eb072ce	Trac #26623: update doctests in cone.py to use new "cones" methods.

[orlitzky]

Author: Michael Orlitzky

[orlitzky]

Branch: u/mjo/ticket/26623

[videlec]

comment:9

I am fine with the class names but the objects would better be like graphs (and many other constructors): that is in caml case like Trivial, NonnegativeOrthant, etc

[orlitzky]

comment:10

Replying to @videlec:

I am fine with the class names but the objects would better be like graphs (and many other constructors): that is in caml case like Trivial, NonnegativeOrthant, etc

Doesn't that violate PEP8? They're not really constructors, even though they do mimic constructors (in the sense that they would be constructors if I didn't want the cones. prefix). However this is python where almost every function creates a new object and returns it, so the distinction seems a little arbitrary to me.

[videlec]

comment:11

Please have a look at: graphs., designs., polytopes., manifolds., ... All of them use caml case and most of them are not classes. I am not talking about PEP8 but consistency within sage.

[novoselt]

comment:12

It seems that when you do specify a lattice, you also have to specify its dimension and there is a check that you do not make a mistake. It would be more pleasant to provide only one input argument: either a number to get the cone in the default lattice of that dimension, or a lattice and get a cone in that lattice.

[orlitzky]

comment:13

Replying to @novoselt:

It seems that when you do specify a lattice, you also have to specify its dimension and there is a check that you do not make a mistake. It would be more pleasant to provide only one input argument: either a number to get the cone in the default lattice of that dimension, or a lattice and get a cone in that lattice.

But that's exactly how Cone() works, and Cone() is the function that underlies all of these convenience methods.

The additional check is not strictly required; it serves only to provide a more relevant error message. For example, if I were to delete the check and then pass a lattice of the wrong rank to Cone(), I would get

sage: M = ToricLattice(2)
sage: Cone([(1,-1,0),(0,1,-1)], lattice=M)
...
TypeError: cannot convert (1, -1, 0) to Vector space of dimension 2 over Rational Field!

If you know the generators of (say) the Schur cone, then you might be able to figure out what that error message means. But, if you don't know what those generators are off the top of your head, then an error message complaining about converting (1,-1,0) to a rational vector space could be confusing. That said, I'm not too attached to the additional sanity check and error message, and would be happy to remove them.

It would be more pleasant to provide only one input argument: either a number to get the cone in the default lattice of that dimension

I don't think this is workable, since in applications we often want to intersect the dual of some cone (which will live in a non-default lattice) with e.g. a nonnegative orthant.

only one input argument... a lattice, and get a cone in that lattice.

This is the best API, but as a user-interface is pretty annoying to use. The "pass a lattice, but only if you don't want the default one" approach of Cone() is IMO a good compromise. Having to do e.g.

sage: cones.trivial(ToricLattice(3))

is even worse than

sage: Cone([], ToricLattice(3))

In the common case, the user is left wondering why he has to type ToricLattice(n) all the time instead of just n. The default is usually fine until you start doing complicated stuff.

[novoselt]

comment:14

Cone works for arbitrary cones and it may be convenient to give input as integer vectors, especially as a user input, in which case it is handy to specify some lattice and cast all vectors to it. Just specifying a lattice does not make much sense for Cone except for a trivial one, perhaps. But in your constructors everything is determined by dimension, there are no explicit rays to give, so what I propose is to have two options: cones.trivial(3) and cones.trivial(my_lattice). I do not suggest that you eliminate the form when only the integer is given, rather do not require this integer in those cases, when the lattice was given explicitly!

[orlitzky]

comment:15

Replying to @novoselt:

I do not suggest that you eliminate the form when only the integer is given, rather do not require this integer in those cases, when the lattice was given explicitly!

Ah, sorry I misunderstood. This sounds like a good idea.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from eb072ce to 5533bd2

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

53f7b42	Trac #26623: add global entries for pending sage.geometry.cones citations.
cab6be8	Trac #26623: add cones. shortcuts for common cones.
f48f779	Trac #26623: add the sage.geometry.cones module to the documentation index.
eda8548	Trac #26623: update doctests in cone.py to use new "cones" methods.
8122dc0	Trac #26623: factor out common "cones" argument processing.
5533bd2	Trac #26623: use Pascal case for "cones" method names.

[orlitzky]

comment:17

Each method will now make an attempt to determine the ambient dimension from the lattice and vice-versa. An error is raised if they are both left unspecified, or if they disagree. I don't think it's necessary to prohibit e.g. nonnegative_orthant(3,M) if the rank of M is 3.

I also added a commit on top of the others to switch to PascalCase names. My mother taught me that "everyone else is doing it" isn't a good excuse, so I'm not sold on the motivation, but it's there. I don't feel that strongly about it if I'm in the minority.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from 5533bd2 to 59fd912

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

7ddd803	Trac #26623: add global entries for pending sage.geometry.cones citations.
0588289	Trac #26623: add cones. shortcuts for common cones.
a80699c	Trac #26623: add the sage.geometry.cones module to the documentation index.
1f8aeaf	Trac #26623: update doctests in cone.py to use new "cones" methods.
d266a33	Trac #26623: factor out common "cones" argument processing.
59fd912	Trac #26623: use Pascal case for "cones" method names.

[orlitzky]

comment:19

Force-pushed a rebase onto the "develop" branch to fix a conflict in src/doc/en/reference/references/index.rst.

[jplab]

comment:20

There seems to be a conflict with version 8.9.beta8.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

4fcb2c4	Trac #26623: add global entries for pending sage.geometry.cones citations.
d163d52	Trac #26623: add cones. shortcuts for common cones.
58e305f	Trac #26623: add the sage.geometry.cones module to the documentation index.
65f70a4	Trac #26623: update doctests in cone.py to use new "cones" methods.
b478929	Trac #26623: factor out common "cones" argument processing.
f4c577e	Trac #26623: use Pascal case for "cones" method names.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from 59fd912 to f4c577e

[orlitzky]

comment:22

Replying to @jplab:

There seems to be a conflict with version 8.9.beta8.

Another conflict in the global references list. Fixed now, thanks for the heads up.

[kliem]

comment:23

Would it make sense to use lazy imports?

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from f4c577e to 23b7dfc

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

146bf23	Trac #26623: add global entries for pending sage.geometry.cones citations.
f8f5feb	Trac #26623: add cones. shortcuts for common cones.
0d2d9a6	Trac #26623: add the sage.geometry.cones module to the documentation index.
29b269d	Trac #26623: update doctests in cone.py to use new "cones" methods.
5665a76	Trac #26623: factor out common "cones" argument processing.
b919f90	Trac #26623: use Pascal case for "cones" method names.
23b7dfc	Trac #26623: use a lazy_import for the global "cones" object.

[orlitzky]

comment:25

Replying to @kliem:

Would it make sense to use lazy imports?

Good idea. I dislike lazy imports when they're used to mis-structure a hierarchy, but in this case it prevents sage from fully loading an esoteric feature and should speed things up.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

2cffd93	Trac #26623: add global entries for pending sage.geometry.cones citations.
d2e22a2	Trac #26623: add cones. shortcuts for common cones.
211477f	Trac #26623: add the sage.geometry.cones module to the documentation index.
e707e63	Trac #26623: update doctests in cone.py to use new "cones" methods.
f0542f2	Trac #26623: factor out common "cones" argument processing.
0f4d979	Trac #26623: use Pascal case for "cones" method names.
9d605a7	Trac #26623: use a lazy_import for the global "cones" object.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from 23b7dfc to 9d605a7

[orlitzky]

comment:27

Rebased onto develop branch.

[jplab]

comment:28

Some minor things:

• Is CamelCase really necessary (compare the polytopes. library...)?
• The function _preprocess_args needs doctests.
• There is a pyflake import error, could you take care of it?

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from 9d605a7 to 3c88088

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Branch pushed to git repo; I updated commit sha1. New commits:

3c88088

Trac #26623: remove two unused imports.

[orlitzky]

comment:30

Replying to @jplab:

Some minor things:

• Is CamelCase really necessary (compare the polytopes. library...)?

The hardest question to answer. I started with e.g. cones.nonnegative_orthant(), and was asked to change it for consistency (read the earlier comments). I don't have strong feelings either way, but slightly prefer the way the lowercase/underscore names look. We should have a consistent, documented policy about this. I started a thread on -dev about it.

• The function _preprocess_args needs doctests.

This is doctested, it just doesn't look like it. That code was originally inlined into each individual method, and each individual method still has tests to make sure that it works. You can check e.g. the last three TESTS for Schur() to see that they're really testing the code in _preprocess_args. This actually results in some duplicated tests, but makes sure that I haven't totally forgotten to call _preprocess_args for some cone (and is nice documentation).

• There is a pyflake import error, could you take care of it?

Fixed, thanks.

[orlitzky]

Changed commit from 3c88088 to f4a4d34

[orlitzky]

Changed branch from u/mjo/ticket/26623 to u/mjo/ticket/26623-ng

[orlitzky]

comment:31

If no one else, I have convinced myself that lowercase-with-underscores is the better option here.

Things got messy when I tried to revert the camel-case commit, so I wound up re-doing everything. The new implementation is slightly different:

• I renamed the module from cones to cone_catalog which fits better with the existing naming scheme for algebras, designs, etc. Someone pointed this out on the mailing list and it makes sense to me.
• I got rid of the class, and put everything back into functions exposed at the top level. The whole module is now imported under the global name cones. The reason for not doing this in the first place is that you don't want cones.<tab> to show things like ZZ and matrix that just happen to be imported into the module. After moving a few imports into the functions that use them, that's no longer an issue here.
• Things are back to lowercase-with-underscores.
• The lazy import isn't needed any longer, because we're not at risk of importing a pre-created object that will waste memory. Now it's just function names that do nothing unless you call them.
• All of the commit messages have been updated to reflect these changes.

[kliem]

comment:32

Some minor comments.

• Your INPUT item all end with a period
• You describe twice (after INPUT and after OUTPUT), when errors are being raised.
• You misspelled rearrangenent at least once.
• Please provide doctests for _preprocess_args.
• nonnegative_orthant:
  • Can you be more precise with the reference?

+    REFERENCES:
+
+    - [BV2009]_

-

• rearrangement:
  • "rearrange," (comma in quotes)
  • again the references are rather unspecific except for [Jeong2017]_ (is there any way to obtain [Jeong2017]_?)
  • I would formulate the statement exactly the other way around. I guess figuring out the Lyapunov rank of the nonnegative orthant or a halfspace is the easy part.

+    Jeong's Corollary 5.2.4 [Jeong2017]_ states that if ``p`` is
+    either one or ``ambient_dim``, then the Lyapunov rank of the
+    rearrangement cone in ``ambient_dim`` dimensions is
+    ``ambient_dim``. Moreover for all other values of ``p``, its
+    Lyapunov rank is one::

• Isn't this clear from the definition anyway?

+    Jeong's Proposition 5.2.1 [Jeong2017]_ states that the rearrangenent
+    cone of order ``p`` is contained in the rearrangement cone of
+    order ``p + 1``::

So maybe just the following would suffice

+    The rearrangement cone of order ``p`` is
+    by definition contained in the one of order ``p + 1``::

-

+    The generators for the rearrangement cone are given by [Jeong2017]_
+    Theorem 5.2.3.

Maybe it's more useful for the reader to specify the generators concretely here.

• schur
  • Please provide a definition.
  • Again the references would be more helpful with concrete information (such as sections etc.)

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Branch pushed to git repo; I updated commit sha1. Last 10 new commits:

ed90795	Trac #26623: reduce INPUT/OUTPUT documentation duplication.
d83f5e0	Trac #26623: add more documentation and examples for _preprocess_args().
06df9d9	Trac #26623: remove an unnecessary citation.
bde9302	Trac #26623: improve the ALGORITHM block for the rearrangement() cone.
6ceefa8	Trac #26623: clarify the example of Jeong's Corollary 5.2.4.
ab5ea7d	Trac #26623: fix a bound in a doctest (off-by-two from range()).
f54d45a	Trac #26623: improve the references for the rearrangement cone.
26618d4	Trac #26623: improve documentation for the schur() cone.
8034754	Trac #26623: remove two breaking spaces in LaTeX code.
42e9e2a	Trac #26623: remove periods from bullet list items.

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from f4a4d34 to 42e9e2a

[orlitzky]

comment:34

Thanks for the feedback. I think I addressed all of your suggestions. I left each change separate for ease of review, but these commits can all be squashed eventually.

The only way I know of to get Juyoung's thesis is to email him and ask. I don't know of another published reference for the Lyapunov rank result, though.

[kliem]

comment:35

Thanks. Overall this will be a nice improvement.

A few comments for now.

• I would suggest using commas here:

+globally-available ``cones`` prefix, to create some common cones:
+
+- The nonnegative orthant,
+- the rearrangement cone of order ``p``,
+- the Schur cone,
+- the trivial cone.
-   globally-available ``cones`` prefix, to create some common cones:
-
-- The nonnegative orthant
-- The rearrangement cone of order ``p``
-- The Schur cone
-- The trivial cone

• The warning at the beginning should include mentioning, that any import will be imported at startup.
• If you lazy import Cone and random_cone than the bots will be fine with the number of startup modules. Also this should be done anyway, shouldn't it?

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

Changed commit from 42e9e2a to 282a694

[7ed8c4ca-6d56-4ae9-953a-41e42b4ed313]

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

fdfd1e4	Trac #26623: add cones. shortcuts for common cones.
ffa1e9b	Trac #26623: add sage.geometry.cone_catalog to the documentation index.
282a694	Trac #26623: update sage.geometry.cone doctests to use the cone catalog.

[orlitzky]

comment:37

Thanks, I added the commas, and expanded the warning at the beginning of the module.

As for the bots: I don't think they're complaining about Cone or random_cone:

1. Both of those are loaded by sage.geometry.all anyway.

2. Those import statements aren't executed until you run one of the functions in the cone_catalog module, so importing them is essentially "lazy" even if you're using sage as a library and not interactively.

Instead, I think the bots are pointing out that sage.geometry.cone_catalog is now loaded into the global namespace as cones, but that's almost unavoidable since that's the whole point. However, it is indeed slower than the lazy import of the class was, should you decide not to use it.

There may be some clever way to delay even that import, while still retaining the name "cones" along with its tab-completion list. For example, if I give a top-level variable name to the module, and then import that object lazily... e.g. this seems to work in cone_catalog.py

import sys

# A python object that can be imported lazily with lazy_import(...).
cones = sys.modules[__name__]

def __dir__():
    # Don't expose any global imports or variables to tab-completion.
    return ['nonnegative_orthant',
            'rearrangement',
            'schur',
            'trivial']

if I then do

lazy_import('sage.geometry.cone_catalog', 'cones')

in sage/geometry/all.py. But that's veering dangerously close to "too clever" to save a few microseconds in my opinion.

---

New commits:

fdfd1e4	Trac #26623: add cones. shortcuts for common cones.
ffa1e9b	Trac #26623: add sage.geometry.cone_catalog to the documentation index.
282a694	Trac #26623: update sage.geometry.cone doctests to use the cone catalog.