Opened 6 years ago
Closed 6 years ago
#18213 closed defect (fixed)
A lot of polytope constructors are broken
Reported by:  vdelecroix  Owned by:  

Priority:  major  Milestone:  sage6.7 
Component:  geometry  Keywords:  
Cc:  ncohen  Merged in:  
Authors:  Vincent Delecroix  Reviewers:  Nathann Cohen 
Report Upstream:  N/A  Work issues:  
Branch:  a58da00 (Commits)  Commit:  a58da00e4907251c2ab47d7af1109da31901ee36 
Dependencies:  #18211  Stopgaps: 
Description (last modified by )
A lot of polytopes constructors in sage.geometry.polyhedron.library
. For example
 The Rhombicuboctahedron does not return what it should (as Python division was thought as Sage integer divisions). There is in the code
verts = [ [3/2, 1/2, 1/2], [3/2, 1/2, 1/2], [3/2, 1/2, 1/2], ...
 The
great_rhombicuboctahedron
is defined overQQ
but it should be defined overQQ[sqrt(2)]
! There are in two places rough approximation of sqrt(2) in the codev1 = QQ(131739771357/54568400000) # ~ sqrt(2) + 1 but actually 2 due to Python division v2 = QQ(104455571357/27284200000) # ~ 2 sqrt(2) but actually 3 due to Python division
Instead, we should use the base_ring
argument (with appropriate defaults) and use base_ring(2).sqrt()
instead.
 The functions unrelated to construction of Polytopes are moved out of the class
Polytopes
:Polytopes.orthonormal_1
,Polytopes.project_1
will be renamed respectivelyzero_sum_projection
andproject_points
 the one line
Polytopes._pfunc
is just removed
While we're at it, remove the deprecations from #11634.
During the process I discovered two annoying bugs:
 #18214 Bug in volume computation of polyhedron
 #18220 Bug when creating a polyhedron with coefficients in RR
Moreover, we can get a great speed up with the following because many polytopes have now coordinates in a quadratic number fields:
Change History (52)
comment:1 in reply to: ↑ description Changed 6 years ago by
 Description modified (diff)
comment:2 Changed 6 years ago by
 Cc ncohen added
comment:3 Changed 6 years ago by
comment:4 Changed 6 years ago by
 Branch set to u/vdelecroix/18123
 Commit set to e598e289689ad8c294500303c89d5c1a0de2a060
 Description modified (diff)
 Status changed from new to needs_review
New commits:
5110cc9  Trac 18211: use LattE to compute the Ehrhart polynomial

cf2a10f  Trac 18211: another example using Birkhoff_polytope

ca8da5a  Trac 18211: fixed spelling LattE integral*e*

d52053c  Trac 18211: check error + options

0f5122b  Trac 18211: fix command line arguments

b877932  Trac 18211: fix arguments + doctests + cwd=SAGE_TMP

3b30bc8  Trac 18213: rewrite the polyhedron library

e598e28  Trac 18123: fix doctests

comment:5 followup: ↓ 6 Changed 6 years ago by
Depends on #18211?
comment:6 in reply to: ↑ 5 Changed 6 years ago by
 Branch changed from u/vdelecroix/18123 to u/vdelecroix/18213
 Commit changed from e598e289689ad8c294500303c89d5c1a0de2a060 to c44b53efbe0fa899951481b96c4ab0a15ebb88e3
 Dependencies set to #18211
comment:7 Changed 6 years ago by
 Description modified (diff)
comment:8 Changed 6 years ago by
 Commit changed from c44b53efbe0fa899951481b96c4ab0a15ebb88e3 to 6e7cf214ea1ee23be11060e0327d20e2ff1a35c8
Branch pushed to git repo; I updated commit sha1. New commits:
6e7cf21  trac #18211: Seealsos

comment:9 Changed 6 years ago by
 Commit changed from 6e7cf214ea1ee23be11060e0327d20e2ff1a35c8 to eef42afa3e3229534f2ae2d3303fbb176eaecdd8
comment:10 Changed 6 years ago by
 Description modified (diff)
comment:11 followup: ↓ 12 Changed 6 years ago by
The golden_ratio()
can obviously be defined in terms of sqrt(5)
, so I don't see the need for that. I'm also really wondering whether you need the custom sqrt()
function, doesn't QuadraticField(n)
work?
comment:12 in reply to: ↑ 11 ; followup: ↓ 13 Changed 6 years ago by
 Status changed from needs_review to needs_work
Replying to jdemeyer:
The
golden_ratio()
can obviously be defined in terms ofsqrt(5)
, so I don't see the need for that. I'm also really wondering whether you need the customsqrt()
function, doesn'tQuadraticField(n)
work?
Right. The problem of embedding should be done elsewhere (my version embedds in QQbar
whereas the generic QuadraticField(n)
embedds into RLF
).
comment:13 in reply to: ↑ 12 Changed 6 years ago by
Replying to vdelecroix:
Replying to jdemeyer:
The
golden_ratio()
can obviously be defined in terms ofsqrt(5)
, so I don't see the need for that. I'm also really wondering whether you need the customsqrt()
function, doesn'tQuadraticField(n)
work?Right. The problem of embedding should be done elsewhere (my version embedds in
QQbar
whereas the genericQuadraticField(n)
embedds intoRLF
).
Here is one reason. The generator name of QuadraticField(2)
is a
and it better be sqrt2
in this case. For the golden ratio, the advantage is that its get printed as phi
. In particular you have
sage: polytopes.icosahedron().volume() 5/6*phi + 5/6
The answer 5/6*a + 5/6
would have been terrible here. I have nothing against using phi = (sqrt5+1)/2
. So I will at least remove this one. For the quadratic fields, one option is to globally make this change in another ticket. What do you think?
comment:14 Changed 6 years ago by
 Commit changed from eef42afa3e3229534f2ae2d3303fbb176eaecdd8 to c1804d5879eedab9140d07096bfd1ad29ae75fa0
Branch pushed to git repo; I updated commit sha1. New commits:
c1804d5  Trac 18213: (review) remove sqrt/golden_ratio + doc

comment:15 Changed 6 years ago by
 Status changed from needs_work to needs_review
All right. I removed the sqrt
and golden_ratio
function. Needs review!
Vincent
comment:16 followup: ↓ 17 Changed 6 years ago by
Vincent: the diff of your main commit is very, very very hard to read. Could you split in order to make it easier to understand?
comment:17 in reply to: ↑ 16 Changed 6 years ago by
Replying to ncohen:
Vincent: the diff of your main commit is very, very very hard to read. Could you split in order to make it easier to understand?
Not sure. It is basically a complete rewrite. I added a ton of documentation, and modify many functions to fit with the actual definition!
comment:18 followup: ↓ 19 Changed 6 years ago by
Okay. Then could you at least change the description of this ticket to explain all changes that you made? (you apparently reated new functions, for example.. Why?)
comment:19 in reply to: ↑ 18 Changed 6 years ago by
 Description modified (diff)
Replying to ncohen:
Okay. Then could you at least change the description of this ticket to explain all changes that you made? (you apparently reated new functions, for example.. Why?)
done in the ticket description.
comment:20 followup: ↓ 21 Changed 6 years ago by
About golden_ratio
and sqrt
: shouldn't they be in number_field_element_quadratic
instead ?
About zero_sum_projection
what about having it in matrix.<tab>
?
Nathann
comment:21 in reply to: ↑ 20 Changed 6 years ago by
Replying to ncohen:
About
golden_ratio
andsqrt
: shouldn't they be innumber_field_element_quadratic
instead ?
Removed in c1804d5.
About
zero_sum_projection
what about having it inmatrix.<tab>
?
It is very specific. I am not sure how I would call it in matrix.<tab>
.
comment:22 Changed 6 years ago by
 Commit changed from c1804d5879eedab9140d07096bfd1ad29ae75fa0 to 524f00ea3d1d3734a9915cd6052d0bf677fdb6e3
Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:
a9ce461  trac #18211: print stderr when verbose=False

9307916  Trac 18211: "in in" > "in"

a6d03cd  Trac 18211: any option to the command (using **kwds)

c935177  Trac 18211: handle bad output from count

7215956  Trac 18211: write explicitely the options in function declaration

74ef0f1  trac #18211: Nothing important

ec0e946  Trac 18213: merge sage6.7.beta1

0b87ff0  Trac 18213: rewrite the polyhedron library

2ddcdef  Trac 18213: fix doctests

524f00e  Trac 18213: Seealsos

comment:23 Changed 6 years ago by
O_o
comment:24 Changed 6 years ago by
Cleaner version based on sage6.7.beta1 and the complete #18211.
comment:25 followup: ↓ 26 Changed 6 years ago by
Here are some comments. Please understand that doing everything at once makes things *MUCH* harder to understand for the reviewer.
Improving the implementation and improving the doc could have been two separate commits. The work that you do on each function could have been a separate commit. Removing deprecations could have been a separate commit. There are also things that anybody can review (only code) and things for which you need to be used to the tools you deal with.
 The projection is not the orthogona projection I expected:
sage: sage: zero_sum_projection(3)*vector([1,0,0]) (0.7071067811865475, 0.4082482904638631)
You can either change the matrix or emphasize in the docstring that the projection is not unique. Something like "returns a (d1)xd matrix of rank d whose kernel contain (1,...,1). At first sight, I expected the function to return a matrix of dimensionsd\times d
. I updated the docstring on this matter.
 I do not understand the paragraph in the doc of
project_points
but that's probably because I don't know the maths behind.
 Compatible?
This projection is the only one compatible with the restriction to the first coordinates
 You should probably add a
# tol
after thisIts volume is `\sqrt{d+1} / d!`
Or even better, check that the difference between the two is 0 (still, with some tolerance)
 At this point I still do not know if it is very wise to use this definition of the projection, but it would be cool if every function which uses this projection (like !Simplex) could redirect toward the matrix function. This way people would have a chance to learn what exactly they get as a result.
 What about this error?
sage: polytopes.icosahedron(base_ring=QQ) ... TypeError: unable to convert 1/4*sqrt(5) + 1/4 to a rational
I added a small commit at public/18213
Nathann
comment:26 in reply to: ↑ 25 ; followup: ↓ 28 Changed 6 years ago by
Replying to ncohen:
Here are some comments. Please understand that doing everything at once makes things *MUCH* harder to understand for the reviewer.
I am very sorry.
 The projection is not the orthogona projection I expected:
Which one did you expected? I improved the documentation.
 At this point I still do not know if it is very wise to use this definition of the projection, but it would be cool if every function which uses this projection (like !Simplex) could redirect toward the matrix function. This way people would have a chance to learn what exactly they get as a result.
I see. But I am really not confortable in moving this to matrix.<tab>
as well. One possibility is to move it back to polytopes.<tab>
. What do you think?
 What about this error?
sage: polytopes.icosahedron(base_ring=QQ) ... TypeError: unable to convert 1/4*sqrt(5) + 1/4 to a rational
Is it not clear enough? The coordinates of the icosahedron need to be defined in QQ[sqrt(5)]
, hence the error. The following is ok (but very slow, ~3secs on my computer)
sage: I = polytopes.icosahedron(base_ring=AA)
Apart adding your example to the documentation, I do not see much what I can do (see the last commit). A deprecation?
Vincent
comment:27 Changed 6 years ago by
 Commit changed from 524f00ea3d1d3734a9915cd6052d0bf677fdb6e3 to fac3cffdb326dad1d3c61045d1c33044c99b8b32
comment:28 in reply to: ↑ 26 ; followup: ↓ 31 Changed 6 years ago by
Hello,
Which one did you expected? I improved the documentation.
Well, I expected an orthogonal projection on this hyperplane. And I probably expected the base of the hyperplane to be the projection of d1 vectors of the base from the original space I guess.
 At this point I still do not know if it is very wise to use this definition of the projection, but it would be cool if every function which uses this projection (like !Simplex) could redirect toward the matrix function. This way people would have a chance to learn what exactly they get as a result.
I see. But I am really not confortable in moving this to
matrix.<tab>
as well. One possibility is to move it back topolytopes.<tab>
. What do you think?
All I was suggesting in the message above was to add a link in the doc. The message you added in project_points
is perfect
The projection used is the matrix given by :func:`zero_sum_projection`.
To me it should appear whenever there is a 'project' argument in the function.
Is it not clear enough? The coordinates of the icosahedron need to be defined in
QQ[sqrt(5)]
, hence the error.
Precisely: Could you add somewhere in the doc that the field must contain QQ(sqrt(5))
?
Nathann
comment:29 followup: ↓ 30 Changed 6 years ago by
the function dodecahedron does nothing with its base_ring
argument. Also, why should the projection function return 'None' instead of '[]' when it is called with no argument?
Nathann
comment:30 in reply to: ↑ 29 ; followup: ↓ 32 Changed 6 years ago by
Replying to ncohen:
the function dodecahedron does nothing with its
base_ring
argument. Also, why should the projection function return 'None' instead of '[]' when it is called with no argument?
You corrected it in 7eb2326
comment:31 in reply to: ↑ 28 Changed 6 years ago by
Replying to ncohen:
Hello,
Which one did you expected? I improved the documentation.
Well, I expected an orthogonal projection on this hyperplane. And I probably expected the base of the hyperplane to be the projection of d1 vectors of the base from the original space I guess.
The description is explicit in the doc (commit fac3cff)
 At this point I still do not know if it is very wise to use this definition of the projection, but it would be cool if every function which uses this projection (like !Simplex) could redirect toward the matrix function. This way people would have a chance to learn what exactly they get as a result.
I see. But I am really not confortable in moving this to
matrix.<tab>
as well. One possibility is to move it back topolytopes.<tab>
. What do you think?All I was suggesting in the message above was to add a link in the doc. The message you added in
project_points
is perfectThe projection used is the matrix given by :func:`zero_sum_projection`.To me it should appear whenever there is a 'project' argument in the function.
Will do!
Is it not clear enough? The coordinates of the icosahedron need to be defined in
QQ[sqrt(5)]
, hence the error.Precisely: Could you add somewhere in the doc that the field must contain
QQ(sqrt(5))
?
There is one example at the end of the doc (from fac3cff)
comment:32 in reply to: ↑ 30 Changed 6 years ago by
You corrected it in 7eb2326
Arggggg... I was looking at the original commit again. I was wondering why it has been removed >_<
The description is explicit in the doc (commit fac3cff)
Yes yes now it's good!
There is one example at the end of the doc (from fac3cff)
It happens several times though. Could you add it in the doc of base_ring
? You already talk about this field there. It could be something like "If set to None then it will be QQ[\phi]
, otherwise it must be a field that contains it".
Nathann
comment:33 Changed 6 years ago by
The 600cell code would look a bit better if you were building all points at once instead of interrupting the flow to decide which ring you should work with.
comment:34 Changed 6 years ago by
Funny. There is no automorphism_group
function for polytopes O_o
comment:35 Changed 6 years ago by
OUTPUT: EXAMPLES::
comment:36 Changed 6 years ago by
Replying to ncohen:
Funny. There is no
automorphism_group
function for polytopesO_o
Yes! This is very sad... But you have to be careful about what it means. There is the combinatorial automorphism group and the isometry group. I have no idea how to implement the latter efficiently.
comment:37 Changed 6 years ago by
Most probably the ugliest thing I ever saw
sage: Sequence([Graph()]).universe() <class 'sage.graphs.graph.Graph'> sage: Sequence([Graph(),1]).universe() Category of objects sage: Sequence([1]).universe() Integer Ring
comment:38 Changed 6 years ago by
HMmmmmmmm... I don't exactly know what the isometry group is, but let's try anyway: what about defining a complete graph on your points, in which each edge has a color associated to its length. Wouldn't the automorphism group of that be what you want?
comment:39 Changed 6 years ago by
Okayyyyyyyyyyyyyy. End of the review. Nothing else to add ;)
Very good job. This code needed to be cleaned, and it is much better now.
But please, next time: smaller patches. You waste your reviewers' health :P
Nathann
comment:40 Changed 6 years ago by
Replying to ncohen:
HMmmmmmmm... I don't exactly know what the isometry group is, but let's try anyway: what about defining a complete graph on your points, in which each edge has a color associated to its length. Wouldn't the automorphism group of that be what you want?
That is a valid definition. But the standard one is the set of orthogonal matrices that preserve the polyhedron. And it would be nice for isometry_group
to be a finite matrix group (and not permutations)!
comment:41 Changed 6 years ago by
 Commit changed from fac3cffdb326dad1d3c61045d1c33044c99b8b32 to 1a82754ea33046384dbefc6045de5bef0e0a2261
comment:42 Changed 6 years ago by
Patchbot should be happy after that.
comment:43 followup: ↓ 48 Changed 6 years ago by
Some broken doctests in the patchbot's report. Also, can you be sure that the sqrt(6.) / factorial(5)
are not subject to numerical noise? There is no #tol
on them right now.
Nathann
comment:44 Changed 6 years ago by
Hmmmmmm.. I can't seem to find on google an implementation of an isometry group function anywhere O_o
comment:45 Changed 6 years ago by
Hi Vincent and Nathann
Just a quick comment, as I'm not sure how useful this will be for general polyhedra: There is a method for lattice polytopes:
LatticePolytope_PPL
(in sage.geometry.polyhedron.ppl_lattice_polytope
) has a method lattice_automorphism_group()
Cheers, Jan
comment:46 Changed 6 years ago by
 Status changed from needs_review to needs_work
comment:47 Changed 6 years ago by
 Commit changed from 1a82754ea33046384dbefc6045de5bef0e0a2261 to a58da00e4907251c2ab47d7af1109da31901ee36
comment:48 in reply to: ↑ 43 Changed 6 years ago by
Replying to ncohen:
Some broken doctests in the patchbot's report. Also, can you be sure that the
sqrt(6.) / factorial(5)
are not subject to numerical noise? There is no#tol
on them right now.
I don't know... Patchbot did not complain about this.
comment:49 Changed 6 years ago by
 Status changed from needs_work to needs_review
comment:50 Changed 6 years ago by
 Reviewers set to Nathann Cohen
 Status changed from needs_review to positive_review
Okayyyyyyy. Then positive_review
!
Nathann
comment:51 Changed 6 years ago by
 Summary changed from A lot of polytopes constructor are broken to A lot of polytope constructors are broken
comment:52 Changed 6 years ago by
 Branch changed from u/vdelecroix/18213 to a58da00e4907251c2ab47d7af1109da31901ee36
 Resolution set to fixed
 Status changed from positive_review to closed
I have a big commit to push in few minutes...