Add the classical construction of the 120-cell

[jplab]

Ticket #27760 introduced the 120-cell in the library, but it does not give the classical construction of Coxeter form 1969 available on Wikipedia:

https://en.wikipedia.org/wiki/120-cell

This ticket provides this classical construction which is much faster since it uses only a degree 2 extension of the rationals.

Depends on #27760

CC: @w-bruns @jplab @mkoeppe @videlec

Component: geometry

Keywords: polytopes, non-rational, normaliz

Author: Jean-Philippe Labbé

Branch/Commit: 9464179

Reviewer: Jonathan Kliem

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

[jplab]

Branch: u/jipilab/classic120

[jplab]

Last 10 new commits:

0c5b8d9	Fix py3 code and tests for Sage8.9beta3
d66f7e2	pep8
a6a5452	fixed some doctests
71fbfdf	Fixed a test
d4260a6	trivial doc fix
95d6bda	fix coding in library
5971228	Added 'polytope' and fixed sum
a1ad959	Made a doctest not tested due to time
cc1a6ab	Merge branch 'develop' into h4_polytopes
14b4918	Added classical construction of 120-cell

[jplab]

Commit: 14b4918

[jplab]

Dependencies: #27760

[sagetrac-git]

Changed commit from 14b4918 to 3ee9569

[sagetrac-git]

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

bd1fd21	fix details
3ee9569	Merge branch 'public/uniform_H4_polytopes' of trac.sagemath.org:sage into h4_polytopes

[jplab]

comment:3

I merged the missing commit from #27760 to get a clean bot test.

[tscrim]

comment:4

I think this can be simplified:

-list(set([tuple(p) for p in flatten([Permutations(list(c)).list() for c in cp])]))
+list(set([tuple(p) for c in cp for p in Permutations(list(c))]))

so there are less transient (large) lists. Also, you should be able to do

verts += itertools.chain.from_iterable([p(tuple(c)) for p in even_perm] for c in cp)

since += for lists works for iterators and you might as well pass a iterable to from_iterable to avoid storing the entire list in memory (twice). Next, I would use latex formatting here: of type `H_4`::. Finally, is it reasonable to do this test with the backend that comes standard with Sage? In that same vein, how long does the normaliz test take?

[jplab]

comment:5

Hi Travis,

Thanks for the feedback, I'll do the changes.

The new implementation with normaliz takes less than 1 sec. Before, it was taking much longer that's why it was tagged not tested.

Now, for the default 'field' backend, I did not test, but I expect it to take a long while. But I'll test it and write it here for the record.

[kliem]

comment:6

A few comments:

• Most importantly exact=False does not work. base_ring is not defined. Also you don't define sqrt5. On that not, it makes sense to define phi outside the if ... else environment.

-            sage: polytopes.one_hundred_twenty_cell(backend='normaliz',construction='as_permutahedron') # not tested - long time
+            sage: polytopes.one_hundred_twenty_cell(backend='normaliz', construction='as_permutahedron') # not tested - long time

# The 64 permutations of [±2,±2,0,0] (the ± are independant)

I count 24.

• I somehow don't like the construction with the Cartesian product. In the comments its all so easy and then I have to think about, why this is exactly true.

  Couldn't we make a use of product from itertools? This creates all possible signs and then we just multiply (elementwise) all permutations:

  verts = [[sign[i]*perm[i] for i in range(4)] for sign in product(*[[1,-1]]*4) for perm in verts_without_sign]
  

For the record:

sage: %time polytopes.one_hundred_twenty_cell(backend='normaliz')
CPU times: user 1.04 s, sys: 8 ms, total: 1.05 s
Wall time: 255 ms
A 4-dimensional polyhedron in (Number Field in sqrt5 with defining polynomial x^2 - 5 with sqrt5 = 2.236067977499790?)^4 defined as the convex hull of 600 vertices
sage: %time polytopes.one_hundred_twenty_cell(backend='field')
CPU times: user 18.1 s, sys: 28 ms, total: 18.2 s
Wall time: 18.2 s
A 4-dimensional polyhedron in (Number Field in sqrt5 with defining polynomial x^2 - 5 with sqrt5 = 2.236067977499790?)^4 defined as the convex hull of 600 vertices
sage: %time polytopes.one_hundred_twenty_cell(backend='normaliz',construction='as_permutahedron')
CPU times: user 17.2 s, sys: 52 ms, total: 17.2 s
Wall time: 15.5 s

[kliem]

comment:7

Replying to @kliem:

• I somehow don't like the construction with the Cartesian product. In the comments its all so easy and then I have to think about, why this is exactly true.

  Couldn't we make a use of product from itertools? This creates all possible signs and then we just multiply (elementwise) all permutations:

  verts = [[sign[i]*perm[i] for i in range(4)] for sign in product(*[[1,-1]]*4) for perm in verts_without_sign]
  

I should have read the code better...
Didn't know that cartesian product does exactly that. Just figured it's way to many lines for something that simple and then I stopped reading and tried to find a shorter solution.

Still one could do all the signs with one line instead of doing it over and over again. Something like:

-            # The 64 permutations of [±2,±2,0,0] (the ± are independant)
-            verts = Permutations([0,0,2,2]).list() + Permutations([0,0,-2,-2]).list() + Permutations([0,0,2,-2]).list()
-
-            # The 64 permutations of the following vectors:
-            # [±1,±1,±1,±sqrt(5)]
-            # [±1/phi^2,±phi,±phi,±phi]
-            # [±1/phi,±1/phi,±1/phi,±phi^2]
-            from sage.categories.cartesian_product import cartesian_product
-            from sage.misc.flatten import flatten
-            full_perm_vectors = [[[1,-1],[1,-1],[1,-1],[-sqrt5,sqrt5]],
-                                 [[phi_inv**2,-phi_inv**2],[phi,-phi],[phi,-phi],[-phi,phi]],
-                                 [[phi_inv,-phi_inv],[phi_inv,-phi_inv],[phi_inv,-phi_inv],[-(phi**2),phi**2]]]
-            for vect in full_perm_vectors:
-                cp = cartesian_product(vect)
-                # The cartesian product creates duplicates, so we reduce it:
-                verts += list(set([tuple(p) for p in flatten([Permutations(list(c)).list() for c in cp])]))
-
-            # The 96 even permutations of [0,±1/phi^2,±1,±phi^2]
-            # The 96 even permutations of [0,±1/phi,±phi,±sqrt(5)]
-            # The 192 even permutations of [±1/phi,±1,±phi,±2]
-            import itertools
-            even_perm_vectors = [[[0],[phi_inv**2,-phi_inv**2],[1,-1],[-(phi**2),phi**2]],
-                                 [[0],[phi_inv,-phi_inv],[phi,-phi],[-sqrt5,sqrt5]],
-                                 [[phi_inv,-phi_inv],[1,-1],[phi,-phi],[-2,2]]]
-            even_perm = AlternatingGroup(4)
-            for vect in even_perm_vectors:
-                cp = cartesian_product(vect)
-                # The cartesian product creates duplicates, so we reduce it:
-                verts += list(itertools.chain.from_iterable([[p(tuple(c)) for p in even_perm] for c in cp]))
+            # The permutations of:
+            # [2,2,0,0]
+            # [1,1,1,sqrt(5)]
+            # [1/phi^2,phi,phi,phi]
+            # [1/phi,1/phi,1/phi,phi^2]
+            vectors = [[2,2,0,0], [1,1,1,sqrt5], [1/phi**2,phi,phi,phi], [1/phi, 1/phi, 1/phi, phi**2]]
+            verts_no_sign = list(x for vec in vectors for x in Permutations(vec).list())
+
+            # The even permutations of [0,1/phi^2,1,phi^2]
+            # The even permutations of [0,1/phi,phi,sqrt(5)]
+            # The even permutations of [1/phi,1,phi,2]
+            import itertools
+            vectors = [[0,1/phi**2,1,phi**2], [0,1/phi,phi,sqrt5], [1/phi,1,phi,2]]
+            even_perm = AlternatingGroup(4)
+            verts_no_sign += list(itertools.chain.from_iterable([[p(tuple(vec)) for p in even_perm] for vec in vectors]))
+
+            # Adding for each vertex copies in all orthants:
+            # [1,1,1,sqrt(5)] -> [±1,±1,±1,±sqrt(5)]
+            from itertools import product
+            verts = [[sign[i]*perm[i] for i in range(4)] for sign in product(*[[1,-1]]*4) for perm in verts_no_sign]

But I guess it's up to taste, what one prefers.

[sagetrac-git]

Changed commit from 3ee9569 to 7defa3c

[sagetrac-git]

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

2f15255	Merge branch 'u/jipilab/classic120' into sage 9.0.beta2
7defa3c	fixed the construction

[jplab]

comment:10

I did the changes. I decided to keep the construction as is.

I added a raise ValueError if one tries to construct it with RDF as cdd returns a numerical inconsistency error.

[kliem]

comment:11

Some minor things:

• The order should be consistent here:

+            # The 64 permutations of [±2,±2,0,0] (the ± are independant)
+            verts = Permutations([0,0,2,2]).list() + Permutations([0,0,-2,-2]).list() + 
Permutations([0,0,2,-2]).list()

also I still count 24 permutations as I already mentioned.

• I'm a bit puzzled, where you remove the duplicates here:

+                # The cartesian product creates duplicates, so we reduce it:
+                verts += itertools.chain.from_iterable([p(tuple(c)) for p in even_perm] for c in cp)

• Instead of importing itertools I would vote for importing just itertools.chain.
  Otherwise you might also consider to use itertools.product instead of sage.categories.cartesian_product.
• Whats the point of naming it K here?

K = QuadraticField(5, 'sqrt5')
base_ring = K

[jplab]

Changed commit from 7defa3c to 5e37f0e

[jplab]

comment:13

Replying to @kliem:

• The order should be consistent here:

+            # The 64 permutations of [±2,±2,0,0] (the ± are independant)
+            verts = Permutations([0,0,2,2]).list() + Permutations([0,0,-2,-2]).list() + 
Permutations([0,0,2,-2]).list()

also I still count 24 permutations as I already mentioned.

Typo. Fixed.

• I'm a bit puzzled, where you remove the duplicates here:

+                # The cartesian product creates duplicates, so we reduce it:
+                verts += itertools.chain.from_iterable([p(tuple(c)) for p in even_perm] for c in cp)

Well, this does what I want, what else can I say? It is very annoying to write something that creates them without duplicates, and we are dealing with a negligible amount of vertices, so I see now reason to make a fuss about inefficiency here.

• Instead of importing itertools I would vote for importing just itertools.chain.

Done.

• Whats the point of naming it K here?

K = QuadraticField(5, 'sqrt5')
base_ring = K

That's a wonderful rhetorical question. My rhetorical answer: to make you ask rhetorical questions. ;-)

---

New commits:

6454b0b	Added classical construction of 120-cell
8b64dec	fixed the construction
5e37f0e	some more fixes

[jplab]

Changed branch from u/jipilab/classic120 to u/jipilab/classic120bis

[kliem]

comment:15

Replying to @jplab:

Replying to @kliem:

• I'm a bit puzzled, where you remove the duplicates here:

+                # The cartesian product creates duplicates, so we reduce it:
+                verts += itertools.chain.from_iterable([p(tuple(c)) for p in even_perm] for c in cp)

Well, this does what I want, what else can I say? It is very annoying to write something that creates them without duplicates, and we are dealing with a negligible amount of vertices, so I see now reason to make a fuss about inefficiency here.

Well, as every entry is unique, the cartesian product creates no duplicates. Hence, you don't have to remove them.

This is why I find this comment very confusing. You comment that you will remove duplicates, but neither are there duplicate nor do you do any instructions that would remove duplicates, if they where there.

Also you have a typo in your chain import, so this doesn't build.

[jplab]

comment:17

Oops, yes, I changed the import to from.

Perhaps my comment is not clear enough, it is not the cartesian product but the action of the even permutations:

+            import itertools import chain
+            even_perm_vectors = [[[0],[phi_inv**2,-phi_inv**2],[1,-1],[-(phi**2),phi**2]],
+                                 [[0],[phi_inv,-phi_inv],[phi,-phi],[-sqrt5,sqrt5]],
+                                 [[phi_inv,-phi_inv],[1,-1],[phi,-phi],[-2,2]]]
+            even_perm = AlternatingGroup(4)
+            for vect in even_perm_vectors:
+                cp = cartesian_product(vect)
+                # The cartesian product creates duplicates, so we reduce it:
+                verts += chain.from_iterable([p(tuple(c)) for p in even_perm] for c in cp)
                                               ^-------p acts here and creates duplicates

[sagetrac-git]

Changed commit from 5e37f0e to 6cb1552

[sagetrac-git]

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

6cb1552

more fixes

[sagetrac-git]

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

9464179

Last edits

[sagetrac-git]

Changed commit from 6cb1552 to 9464179

[jplab]

comment:20

Should be good to go now.

[kliem]

Reviewer: Jonathan Kliem

[kliem]

comment:21

LGTM

[vbraun]

Changed branch from u/jipilab/classic120bis to 9464179