Add Vinberg's algorithm #4023
Add Vinberg's algorithm #4023
Conversation
Codecov ReportAttention: Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## master #4023 +/- ##
==========================================
+ Coverage 84.10% 84.66% +0.55%
==========================================
Files 612 613 +1
Lines 83067 83393 +326
==========================================
+ Hits 69865 70606 +741
+ Misses 13202 12787 -415
|
| series = {Tata Inst. Fundam. Res. Stud. Math.}, | ||
| volume = {No. 7}, | ||
| pages = {323--348}, | ||
| publisher = {Published for the Tata Institute of Fundamental Research, Bombay by Oxford University Press, Bombay}, |
There was a problem hiding this comment.
Hmm, did you run bibtool as described here?
There was a problem hiding this comment.
Also, if there is a DOI available, please add it as well.
| author = {Turkalj, I.}, | ||
| title = {Reflective Lorentzian lattices of signature (5, 1)}, | ||
| year = {2017}, | ||
| school = {TU Dortmund} |
There was a problem hiding this comment.
There is also a paper, could that be cited instead? See https://doi.org/10.1016/j.jalgebra.2018.06.013
If not, a URL or better a DOI to the thesis would be great
src/NumberTheory/vinberg.jl
Outdated
| end | ||
|
|
||
| @doc raw""" | ||
| vinberg_algorithm(Q::ZZMatrix, upper_bound::ZZRingElem; v0::ZZMatrix, root_lengths::Vector{ZZRingElem}, direction_vector::ZZMatrix) -> Vector{ZZMatrix} |
There was a problem hiding this comment.
Q is a matrix, not a lattice, though the text below refers to it as "a lattice". This always opens up ambiguities as to how to interpret that matrix. I thought we had types specifically to represent lattices, such as ZZLat, why not use that?
There was a problem hiding this comment.
Indeed, I agree with fingolfin please prepare a dispatch via ZZLat.
There was a problem hiding this comment.
@ElMerkel You marked this as resolved, but I cannot find the method with ZZLat as argument. Did you upload your changes?
There was a problem hiding this comment.
Oh I did not know that you can see when I mark a conversation as resolved. I thought that I can use this for having a better overall view of the feedback. I am sorry, I can imagine that this led to confusion. The code should be uploaded now.
|
Needs sorting and possibly moved to experimental? @simonbrandhorst any idea? It's your baby....(indirectly) |
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
Co-authored-by: Max Horn <[email protected]>
| if _crystallographic_condition(Q, v) | ||
| if _check_coorientation(Q, roots, v, v0) |
There was a problem hiding this comment.
This is the hot loop. So these two function will be called many times. Can you make them non-allocating?
|
The math looks good to me. |
Co-authored-by: Simon Brandhorst <[email protected]>
Co-authored-by: Simon Brandhorst <[email protected]>
Co-authored-by: Simon Brandhorst <[email protected]>
Co-authored-by: Simon Brandhorst <[email protected]>
Co-authored-by: Simon Brandhorst <[email protected]>
| gcd = reduce(gcd, v) | ||
| if isone(gcd) | ||
| return v | ||
| else | ||
| return (1 // gcd * v) | ||
| end |
There was a problem hiding this comment.
Minor comment: I would recommend not re-using a name like gcd for a local variable, perhaps call it g.
More importantlu, though, I don't think the else part works, or at least it doesn't when I tried it locally:
julia> v = matrix(ZZ, 6 * [1 2 ; 3 4])
[ 6 12]
[18 24]
julia> g = reduce(gcd, v)
6
julia> (1 // g * v)
ERROR: Cannot promote to common type
Even if it works, I would expect this to return a QQMatrix, which is probably not what you want. However, this works:
julia> v / g
[1 2]
[3 4]
So overall perhaps use this?
| gcd = reduce(gcd, v) | |
| if isone(gcd) | |
| return v | |
| else | |
| return (1 // gcd * v) | |
| end | |
| g = reduce(gcd, v) | |
| if isone(g) | |
| return v | |
| else | |
| return v / g | |
| end |
or even
| gcd = reduce(gcd, v) | |
| if isone(gcd) | |
| return v | |
| else | |
| return (1 // gcd * v) | |
| end | |
| g = reduce(gcd, v) | |
| return isone(g) ? v : v / g |
| negative_eigenvalue = x | ||
| break | ||
| end | ||
| stopping_condition =+ |
There was a problem hiding this comment.
Looks as if something is wrong here? Perhaps this was meant?
| stopping_condition =+ | |
| stopping_condition += 1 |
| Eigenvalues = eigenvalues(QQ, Q) | ||
| stopping_condition = 0 | ||
| negative_eigenvalue = 0 | ||
| for x in Eigenvalues | ||
| if x < 0 | ||
| negative_eigenvalue = x | ||
| break | ||
| end | ||
| stopping_condition =+ | ||
| end | ||
| @req stopping_condition < length(Eigenvalues) "Q is not of signature (1, n)" |
There was a problem hiding this comment.
Actually, if I understand this code right, perhaps you could write it like this?
| Eigenvalues = eigenvalues(QQ, Q) | |
| stopping_condition = 0 | |
| negative_eigenvalue = 0 | |
| for x in Eigenvalues | |
| if x < 0 | |
| negative_eigenvalue = x | |
| break | |
| end | |
| stopping_condition =+ | |
| end | |
| @req stopping_condition < length(Eigenvalues) "Q is not of signature (1, n)" | |
| ev = eigenvalues(QQ, Q) | |
| pos = findfirst(is_negative, ev) | |
| @req pos === nothing "Q is not of signature (1, n)" | |
| x = ev[pos] |
I added an implementation of Vinberg's algorithm for computing fundamental roots of Lorentzian lattices.
@simonbrandhorst
@StevellM